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Abstract 

An  emerging  use  of  GPS  is  to  provide  accurate  navigation  information  for 
satellites  in  orbit.  The  GPS  satellites  are  designed  to  provide  service  to  terrestrial  users, 
so  the  antenna  array  points  directly  towards  the  Earth  and  uses  a  narrow  primary 
beamwidth.  Because  GEO  altitudes  are  well  above  the  GPS  constellation,  the  Earth 
occludes  most  of  the  GPS  signals  to  the  satellite.  Decreased  satellite  visibility  is 
debilitating,  as  GPS  navigation  requires  at  least  four  visible  satellites  to  determine 
position.  To  assist  with  the  visibility  problem,  the  receiver  can  look  at  the  GPS  satellite 
transmit  antenna  side  lobes,  but  this  does  not  entirely  solve  the  navigation  problem.  GPS 
measurements  are  inherently  bound  by  receiver  clock  errors.  The  clock  error  must  be 
known  or  estimated  in  order  to  obtain  meaningful  ranging  information.  To  obtain  three- 
dimensional  positioning,  at  least  four  satellites  must  be  tracked  to  solve  for  three 
dimensions  of  position  plus  the  receiver  clock  error. 

A  new  method  for  improving  geostationary  navigation  accuracy  using  GPS  is  to 
correct  the  time  error  by  including  Two-Way  Time  Transfer  (TWTT)  measurements. 
TWTT  is  a  technique  in  which  signals  are  simultaneously  exchanged  between  two  clocks, 
and  is  one  of  the  most  accurate  methods  of  comparing  clocks.  By  effectively  removing 
the  clock  error  between  the  GPS  satellite  and  the  GPS  receiver,  TWTT  allows  meaningful 
information  to  be  gathered  when  less  than  four  GPS  satellites  are  available.  The  results 
show  a  21-38%  improvement  in  the  3-D  RMS  position  accuracy  while  using  TWTT 
between  the  GEO  satellite  and  an  atomic  clock  on  the  ground.  There  was  a  60-70% 
improvement  when  the  clock  on  the  ground  was  synchronized  to  GPS  time. 
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USE  OF  TWO-WAY  TIME  TRANSFER  MEASUREMENTS  TO  IMPROVE 


GEOSTATIONARY  SATELLITE  NAVIGATION 

1  Introduction 


1.1  Overview 

The  Global  Positioning  System  (GPS)  was  developed  by  the  United  States 
Department  of  Defense  to  provide  precise  position,  velocity,  and  time  measurements  to 
military  users  around  the  globe.  Such  information  would  provide  U.S.  military  forces  an 
unparalleled  advantage  among  adversaries,  allowing  accurate  navigation  in  any  weather 
conditions  during  day  or  night  operations.  Upon  initial  operational  capability  on 
December  8,  1991  [17],  civilian  use  of  the  GPS  signal  was  intentionally  corrupted  with 
random  error,  allowing  only  the  military  users  to  have  precise  navigation.  On  May  2, 
2000,  the  full  accuracy  of  civilian  GPS  was  made  freely  available  to  users  of  every  nation 
[16].  Its  role  in  daily  life  has  broadened  at  a  rapid  rate,  and  now  much  of  modern 
infrastructure  is  dependent  upon  GPS.  Because  of  the  growing  reliance  on  GPS,  users 
demand  increased  precision  from  the  system. 

An  emerging  use  of  GPS  is  to  provide  accurate  navigation  information  for 
satellites  in  orbit.  Previously,  radar  and  optical  range  measurements  to  the  satellite  from 
the  ground  were  collected,  allowing  operators  to  determine  the  position  of  the  satellite. 
With  GPS,  the  satellite  can  carry  a  receiver  and  calculate  its  own  position  using  the  GPS 
signals,  eliminating  the  numerous  ground  stations  and  support  costs.  Additionally,  GPS 
navigation  allows  the  satellite  to  operate  with  a  certain  level  of  autonomy  in  station¬ 
keeping,  making  orbital  corrections  without  the  need  for  ground  interaction.  GPS 
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navigation  is  effective  for  satellites  in  Low  Earth  Orbit  (LEO),  but  satellites  operating  in 
Geosynchronous  Orbit  (GEO)  or  Highly  Elliptical  Orbit  (HEO)  experience  extreme 
performance  degradation,  for  reasons  described  below. 

GPS  measurements  are  inherently  bound  by  receiver  clock  errors.  The  clock  error 
must  be  known  or  estimated  in  order  to  obtain  meaningful  ranging  information.  To 
obtain  three-dimensional  positioning,  at  least  four  satellites  must  be  tracked  to  solve  for 
three  dimensions  of  position  plus  the  receiver  clock  error.  Because  the  receiver  clock 
error  can  be  estimated  with  four  measurements,  a  user  can  have  a  very  poor  local  clock 
attached  to  their  receiver  and  still  have  accurate  positioning. 

The  GPS  satellites  (approximately  20,200  kilometers  altitude)  are  designed  to 
provide  service  to  terrestrial  users,  so  the  antenna  array  points  directly  towards  the  Earth 
and  uses  a  narrow  primary  beamwidth.  Because  GEO  (approximately  35,000  kilometers 
altitude)  and  HEO  altitudes  are  well  above  the  GPS  constellation,  the  Earth  occludes 
most  of  the  GPS  signals  to  the  satellite.  Decreased  satellite  visibility  can  be  debilitating, 
as  GPS  navigation  requires  at  least  four  visible  satellites  to  determine  position. 

As  a  solution  to  the  visibility  problem,  the  receiver  can  look  at  more  than  just  the 
primary  GPS  transmit  beam  that  is  intended  for  terrestrial  use.  The  GPS  transmit 
antennas  have  side  lobes,  shown  in  Figure  1.1,  though  at  much  lower  power  levels  than 
the  primary  beam.  These  side  lobes  provide  increased  visibility  for  GEO  and  HEO 
satellites  and  are  a  requirement  for  high-altitude  navigation. 
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GPS  side  lobe 


GPS  side  lobe 


Figure  1.1:  GPS  Signal  Reception  at  GEO  and  HEO  Orbital  Altitudes  [1] 


A  current  method  for  GPS  navigation  in  high-altitude  orbits  is  to  use  high- 
sensitivity  GPS  receivers  and  more  powerful  signal  acquisition  engines  [1],  By 
increasing  the  sensitivity  of  the  receiver,  the  very  weak  GPS  signals  become  usable  and 
the  number  of  visible  satellites  increases,  allowing  for  more  precise  positioning. 

A  new  method  for  improving  geostationary  navigation  accuracy  using  GPS  is  to 
correct  the  time  error  by  including  Two-Way  Time  Transfer  measurements.  Two-Way 
Time  Transfer  (TWTT)  is  a  technique  in  which  signals  are  simultaneously  exchanged 
between  two  clocks.  If  the  distance  between  the  two  clocks  is  equal,  such  as  static 
clocks,  the  propagation  delays  cancel  and  the  difference  between  the  two  clocks  can  be 
precisely  measured  [13].  Dynamic  TWTT  allows  simultaneous  signal  exchanges  to  occur 


3 


between  moving  clocks,  allowing  vehicles,  such  as  airplanes  or  satellites,  to  take 
advantage  of  TWTT  [2],  TWTT  is  one  of  the  most  accurate  methods  of  comparing 
clocks,  so  it  can  be  used  to  create  more  accurate  positioning  using  GPS  by  effectively 
removing  the  clock  error  between  the  GPS  satellite  and  the  GPS  receiver.  Using  TWTT 
measurements,  the  clock  error  could  be  resolved,  allowing  less  than  four  visible  GPS 
satellites  to  provide  meaningful  ranging  information. 

1.2  Related  Research 

1.2.1  GPS  Measurements  Collected  from  Geosynchronous  Transfer  Orbit 

A  group  of  students  at  the  United  States  Air  Force  Academy  designed  and 
constructed  the  “Falcon  Gold”  experiment  that  flew  in  1997.  The  goal  of  this  project  was 
to  measure  GPS  signals  at  high  orbital  altitudes  using  low-cost,  low-power,  off-the-shelf 
components  [19],  The  purpose  was  to  test  concepts,  not  to  collect  high-quality  data. 

The  Falcon  Gold  experiment  was  a  secondary  payload  on  a  DSCS  spacecraft 
launch  and  was  attached  the  side  of  the  Centaur  upper  stage.  Once  the  DSCS  spacecraft 
separated  from  the  upper  stage  and  maneuvered  into  a  geosynchronous  orbit,  the  upper 
stage  continued  to  orbit  in  a  geosynchronous  transfer  orbit  (GTO).  The  GTO  orbit 
apogee  and  perigee  were  35,200  and  200  kilometers,  respectively,  allowing  the  Falcon 
Gold  to  collect  measurement  data  above  and  below  the  GPS  constellation. 

The  experiment  used  an  inexpensive  two-inch  patch  antenna  for  the  GPS  receiver, 
which  was  modeled  and  used  in  this  research  simulation.  Because  of  the  hardware  setup 
and  limited  ground  link  ability,  only  intermittent  data  could  be  collected.  A  total  of  12 
data  frames  were  collected,  detecting  25  PRN  signals,  as  seen  in  Figure  1.2. 
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Figure  1.2:  Falcon  Gold  GPS  Signal  Detection  Locations  [19] 

The  Falcon  Gold  experiment  is  one  of  the  first  examples  of  GPS  signal  detection 
at  altitudes  above  the  GPS  constellation  [19].  The  future  quality  of  performance  can  only 
improve,  since  Falcon  Gold  was  able  to  achieve  these  capabilities  by  using  only  low-cost, 
low-power  hardware.  This  demonstration  was  an  important  step  in  the  progression  of 
GPS-enabled  satellites  in  high-altitude  orbits. 

The  Falcon  Gold  experiment  was  also  able  to  detect  GPS  sidelobe  signals.  This  is 
very  important  for  applications  that  could  benefit  from  the  extended  coverage  offered  by 
sidelobe  signals,  such  as  high-altitude  GPS  positioning.  New  concepts  may  now  be  fully 
explored  without  hesitation,  since  the  physical  act  of  detecting  sidelobes  signals  has  been 
proven  through  experimentation. 

1.2.2  High-Sensitivity  GPS  Receivers  Built  for  GEO  Operation 

The  National  Aeronautics  and  Space  Administration  (NASA)  developed  a  highly 
sensitive  GPS  receiver  for  use  in  geosynchronous  orbits  that  is  designed  to  track  the  very 
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weak  GPS  side  lobe  signals  [1].  This  new  receiver,  called  the  Navigator,  will  drastically 
reduce  GPS  signal  outages,  providing  more  accurate  positioning.  The  Navigator  is  able 
to  track  weak  signals  that  are  not  usable  by  a  standard  GPS  receiver,  as  seen  in  Figure 
1.3.  The  lower  dashed  line  represents  the  sensitivity  threshold  of  the  Navigator,  while  the 
upper  dashed  line  represents  a  standard  space  receiver. 
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Elapsed  time-  since  1 0^5/2001  00:00:00  (hours) 

Figure  1.3:  Simulated  Received  Power  at  GEO  Orbital  Altitude  [1] 

The  Navigator  is  designed  to  provide  autonomous  navigation  for  spacecraft  and 
includes  the  GPS  Enhanced  Onboard  Navigation  Systems  (GEONS)  software.  GEONS 
is  a  powerful  extended  Kalman  filter  that  uses  an  internal  orbital  dynamics  model  in 
conjunction  with  incoming  measurements  to  generate  a  position  solution,  even  when  less 
than  four  GPS  satellites  are  visible  [1].  The  GEONS  orbital  dynamics  model  is 
extremely  sophisticated  and  can  emulate  a  complex  Earth  gravity  field  and  the  forces  due 
to  solar  and  lunar  gravity,  atmospheric  drag,  and  solar-radiation  pressure.  The  Navigator 
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receiver  feeds  measurements  to  the  GEONS,  which  is  able  to  calculate  precise  GEO 
orbits  to  within  10  meters  3D  root  mean  square  as  seen  in  Table  1.1. 


Table  1.1:  GEONS  State  Estimation  Error  Statistics  [1] 

Radial  Intrack  Crosstrack  3D  r.m.s. 


Geostationary  Orbit 

Position  (mi 

-4.424  +  5,640 

-0.546  ±2.265 

0.675  +  2.462 

7.959 

Velocity  (cm/s] 

-0.038  ±  0.076 

-0.133  +  0.051 

-0.001  +  0.025 

0.096 

Low  Earth  Orbit 

Position  (mi 

0.036+  0.747 

0.215+0.995 

0.103  +  0.360 

1.532 

Velocity  (cm/s) 

-0.022  ±  1.613 

0.002  +  0.079 

0.0243+  0.107 

0.210 

The  research 

conducted  in  this 

thesis  is  similar 

to  the  concept  of 

the  Navigator, 

which  supplements  pseudorange  measurements  with  additional  data  to  produce  a  more 
accurate  navigation  solution.  Instead  of  a  high-fidelity  orbital  mechanics  model,  this 
research  uses  TWTT  to  improve  the  GPS  accuracy  at  high  altitudes.  The  Navigator 
requires  specialized  hardware  and  the  GEONS  software,  while  TWTT  measurements 
could  be  used  with  a  very  simple  GPS  receiver  and  low-quality  onboard  clock. 

The  Navigator  simulation  uses  slightly  different  values  for  the  noise  density  and 
the  gain  patterns  of  the  associated  antennas,  so  a  direct  “apples-to-apples”  comparison  is 
not  possible  for  this  thesis  research.  However,  the  close  similarity  offers  a  good 
benchmark  for  validating  the  simulation  used  in  this  thesis. 

1.2.3  New  GPS  Satellite  Architecture  for  High-Altitude  Spacecraft 

The  Aerospace  Corporation  conducted  a  study  in  1998  to  examine  the  effects  of 
various  GPS  satellite  architectures  and  user  processing  options  on  the  accuracy  of  high 
altitude  orbit  determination  using  GPS  [20].  The  study  was  divided  into  two  main 
investigations  to  examine  the  tradeoff  between  GPS  system  modifications  and  user 
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equipment  complexity.  It  includes  several  hypothetical  GPS  satellite  antenna 
configurations  and  receiver  clock  choices. 

Five  antenna  architectures  were  used  in  Monte  Carlo  simulations:  the  Block  II-A 
main  beam  (42°),  the  Block  II-R  main  beam  (38°),  the  Block  II-R  UHF  crosslink  antenna 
(120°),  and  a  hemispherical  back-side  antenna  combined  with  either  the  II-A  or  II-R  main 
beam.  The  three  receiver  clock  choices  were  a  standard  crystal,  an  improved  crystal,  and 
a  standard  atomic  clock. 


Figure  1.4:  Block  II-R  UHF  Crosslink  Antenna  [20] 

The  simulations  did  show  that  GPS  modifications  would  improve  positioning 
accuracy  at  GEO  by  increasing  signal  availability,  but  would  require  a  significant 
increase  in  cost  and  complexity.  Therefore,  a  complete  redesign  of  the  GPS  satellites  to 
assist  GEO  navigation  is  highly  improbable. 
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Figure  1.5:  Block  II-R  with  Hemispherical  Back-side  Antenna  [20] 


The  most  important  result  from  this  study  concerns  the  choice  of  receiver  clock. 
By  using  a  clock  with  greater  accuracy,  the  positioning  error  decreases  significantly  as 
seen  in  Table  1.2.  Such  a  concept  is  the  crux  of  this  thesis  research,  purporting  that 
introducing  a  more  accurate  clock  via  TWTT  can  significantly  improve  GEO  positioning 
using  GPS. 


Table  1.2:  RMS  Position  Errors  -  GPS  Broadcast  Antenna  Configuration  vs  User  Clock  Quality  [20] 


II-R  Main  (38°) 
Only 

II-A  Main  (42°) 
Only 

UHF  Crosslink 
(120°) 

II-R  Main  (38°) 

+  Backside 

II-A  Main  (42°) 

+  Backside 

3D 

16.9  m 

3D 

15.3  m 

3D 

3.6  m 

3D 

4.5  m 

3D 

4.5  m 

Standard 

R 

6.3  m 

R 

6.0  m 

R 

2.4  m 

R 

2.7  m 

R 

2.7  m 

Crystal 

I 

15.4  m 

I 

13.9  m 

I 

2.5  m 

I 

3.4  m 

I 

3.5  m 

C 

3.3  m 

C 

2.1  m 

C 

1.0  m 

C 

1.0  m 

C 

1.0  m 

3D 

11.8m 

3D 

11.0m 

3D 

3.5  m 

3D 

4.4  m 

3D 

4.4  m 

Improved 

R 

4.6  m 

R 

4.4  m 

R 

2.3  m 

R 

2.7  m 

R 

2.6  m 

Crystal 

I 

10.7  m 

I 

9.8  m 

I 

2.5  m 

I 

3.4  m 

I 

3.5  m 

C 

2.0  m 

C 

2.1  m 

C 

1.0  m 

C 

1.0m 

C 

0.9  m 

3D 

8.4  m 

3D 

8.1  m 

3D 

2.8  m 

3D 

3.4  m 

3D 

3.4  m 

Atomic 

R 

2.4  m 

R 

2.3  m 

R 

1.4  m 

R 

1.7m 

R 

1.7  m 

Clock 

I 

7.8  m 

I 

7.6  m 

I 

2.1  m 

I 

2.8  m 

I 

2.8  m 

C 

1.7  m 

C 

1.6  m 

C 

1.0  m 

C 

0.9  m 

C 

0.9  m 

R  =  radial,  I  =  in-track,  C  =  cross-track 
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1.2.4  Two-Way  Time  Transfer  Measurements  used  in  Relative  GPS  Positioning 


A  previous  AFIT  masters  degree  student  researched  the  effects  of  TWTT  on 
networked  differential  GPS  positioning  [8],  Differential  GPS  techniques  provide  relative 
positioning  between  a  mobile  GPS  receiver  and  a  fixed  receiver  with  well-known 
position  information,  as  shown  in  Figure  1.6.  The  fixed  receiver  compares  its  GPS- 
calculated  position  with  its  known  geographic  position  and  determines  the  corrections 
necessary  for  eliminating  the  difference  between  the  two.  These  corrections  can  then  be 
sent  to  a  mobile  receiver,  providing  the  same  amount  of  positioning  accuracy  as  the 
known  location.  In  a  networked  system,  multiple  receivers  exist  in  various  locations, 
such  as  a  swarm  of  Unmanned  Aerial  Vehicles  (UAVs),  and  are  linked  to  the  same  fixed 
reference  station.  Differential  GPS  accuracy  is  on  the  order  of  centimeters  to  meters, 
depending  upon  the  method  used  [21], 


Reference  Station 
a  Kgu-ati  Location 


Figure  1.6:  Real-Time  Differential  GPS  Diagram  [10] 

This  previous  research,  found  in  [8],  used  a  MATLAB®  simulation  to  explore  the 
effects  of  introducing  TWTT  measurements  to  a  network  of  six  mobile  receivers.  The 
TWTT  measurement  accuracy  standard  deviation  varied  between  3,  0.3,  0.03,  and  0.003 
meters.  Smaller  TWTT  measurement  error  results  in  greater  position  measurement 
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accuracy  as  seen  in  Table  1.3,  which  was  created  using  the  data  produced  in  [8].  Clearly, 
there  is  a  significant  benefit  to  using  TWTT  measurements  with  GPS  positioning. 


Table  1.3:  Position  Accuracy  Improvement  When  Using  TWTT 


Scenario 

3-D  Position  RMS 
Value  (m) 

Percent  Increase 
in  Accuracy 

P 

1.587 

p+TWTT(3m) 

1.415 

10.84% 

/)+TWTT(0.3m) 

0.935 

41.08% 

p+TWTT(0.03m) 

0.890 

43.92% 

p+TWTT(0.003m) 

0.889 

43.98% 

p  =  pseudorange  measurement 
TWTT(x)  =  TWTT  measurement  where  x  =  accuracy  level 


The  methodology  found  in  [8]  is  similar  to  the  research  conducted  for  this  thesis. 
The  main  difference  is  that  this  research  examines  TWTT  measurements  used  with  GPS 
positioning  of  geostationary  satellites,  versus  a  network  of  receivers  on  Earth.  Both 
simulations  use  similar  assumptions,  such  as  ignoring  relativistic  effects  and  other  types 
of  error  that  can  be  modeled,  calculated,  and  removed  from  physical  measurements.  The 
glaring  similarity  is  the  idea  that  TWTT  can  drastically  improve  GPS  positioning 
accuracy. 


1.3  Problem  Statement 

The  main  objective  of  this  research  is  to  examine  the  impact  of  adding  TWTT 
measurements  to  geostationary  satellite  positioning  using  GPS  measurements.  By 
constraining  relative  clock  errors,  TWTT  measurements  can  improve  the  positioning 
accuracy. 

By  using  TWTT  measurements,  a  GPS  receiver  in  geostationary  orbit  could 
possibly  calculate  a  good  navigation  solution  without  needing  a  high-sensitivity  receiver. 
If  the  satellite  did  have  a  high-sensitivity  receiver,  the  navigation  solution  would  be  even 
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better.  TWTT  measurements  would  eliminate  the  need  for  a  precise  clock  on  the  satellite 
because  it  would  only  need  a  precise  reference  clock  on  the  ground.  Essentially,  using 
TWTT  with  a  highly  accurate  clock  on  the  ground  and  a  low  quality  clock  on  the  satellite 
would  be  comparable  to  putting  a  highly  accurate  clock  on  the  satellite  itself.  A  highly 
accurate  clock  enables  the  system  to  meaningfully  use  individual  pseudorange 
measurements,  because  the  clock  error  can  be  accurately  modeled  over  time. 

This  research  involves  five  trade  studies  that  quantify  the  benefits  of  using  TWTT 
measurements  with  GPS.  The  first  trade  study  involves  using  ephemeris  data  from 
several  days  to  ensure  that  the  results  are  not  dependent  upon  the  day  that  is  used.  The 
second  trade  study  examines  the  differences  between  a  simple  GPS  receiver  model  and  a 
complex  GPS  receiver  model  to  determine  the  effect  of  simulation  model  accuracy.  The 
simple  GPS  receiver  model  uses  a  simple  cut-off  angle  between  the  satellites  and  the 
Earth  and  a  constant  pseudorange  measurement  noise  value.  The  complex  GPS  receiver 
model  uses  gain  patterns  for  the  GPS  and  geostationary  satellite  antennas  and  a  varying 
pseudorange  measurement  noise  that  depends  on  received  signal  strength.  The  third  trade 
study  examines  the  differences  between  a  rudimentary  gain  pattern  and  a  more  accurate 
gain  pattern  for  the  GPS  satellite  antenna  to  determine  the  effect  of  gain  pattern  accuracy 
levels.  The  fourth  trade  study  investigates  the  effects  of  using  different  types  of  clocks 
for  the  geostationary  satellite  clock  and  the  TWTT  reference  clock  to  determine  the  effect 
of  using  a  high-accuracy  clock  versus  a  low-accuracy  clock.  The  fifth  trade  study 
evaluates  the  results  as  a  function  of  TWTT  measurement  accuracy. 
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1.4  Methodology 

A  simulation  using  MATLAB®  examines  the  potential  benefits  of  using  TWTT 
measurements  in  addition  to  GPS  measurements  for  determining  the  position  of  a 
geostationary  satellite  in  orbit.  The  simulation  is  comprised  of  five  main  functions, 
shown  in  Figure  1.7,  which  involve  loading  parameters,  generating  truth  data,  generating 
measurement  data,  executing  a  Kalman  filter,  and  analyzing  the  results. 


r - 

toad_p 

L _ 

- ^ 

arams 

_ ^ 

r 

generate_truth  * 

l _ j 

r 

generate_meas 

‘v _ / 

r 

*  kalman 

v _ J 

r - * - ^ 

analyze_results 

L ^ 

Figure  1.7:  Simulation  Block  Diagram 


The  load _params  function  allows  the  user  to  input  the  desired  parameters  into  the 
simulation.  These  parameters  are  used  by  the  generate _truth,  generate jneas,  and 
kalman  functions.  The  generate _truth  function  generates  the  ‘true’  values  for  the 
positions  of  each  GPS  satellite  and  the  position  of  the  geostationary  satellite.  The 
generatejneas  function  uses  the  ‘true’  data  to  calculate  pseudorange  measurements  and 
the  TWTT  measurements,  simulating  what  an  actual  GPS  receiver  would  collect.  The 
kalman  function  uses  the  generated  measurements  in  a  Kalman  filter  to  estimate  the 
position  of  the  geostationary  satellite.  The  analyze  Jesuits  function  compares  the  ‘true’ 
data  with  the  Kalman  filter  results  to  determine  how  accurately  the  filter  was  able  to 
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position  the  geostationary  satellite  using  the  measurements  provide  by  generatejneas.  A 
more  detailed  description  of  each  function  is  located  in  Chapter  Three. 


1.5  Thesis  Overview 

Chapter  Two  describes  the  fundamental  topics  that  are  related  to  this  research. 
This  chapter  includes  backgrounds  of  GPS,  orbital  reference  frames,  Two-Way  Time 
Transfer,  and  the  Kalman  fdter.  The  topic  of  GPS  includes  an  overview  of  pseudorange 
measurements  and  associated  equations,  while  the  topic  of  TWTT  includes  an  overview 
of  static  and  dynamic  methods  and  associated  equations.  The  topic  of  orbital  reference 
frames  describes  the  Earth-Centered  Earth-Fixed  (ECEF)  and  Earth-Centered  Inertial 
(ECI)  coordinate  systems. 

Chapter  Three  describes  a  simulation  of  positioning  a  geostationary  satellite  in 
orbit  using  simulated  GPS  measurements  with  and  without  TWTT  measurements.  Using 
the  simulation,  five  trade  studies  were  conducted  to  examine  the  impact  on  system 
performance.  The  trade  studies  include  comparing  the  results  obtained  from  using 
ephemeris  data  from  different  days,  the  differences  between  a  simple  and  complex  GPS 
receiver  model,  the  differences  between  simple  and  complex  gain  pattern  data,  the 
differences  between  using  high-  and  low-accuracy  clocks  for  the  geostationary  and 
TWTT  reference  clocks,  and  the  differences  in  adjusting  the  TWTT  measurement  noise 
error. 

Chapter  Four  discusses  the  results  of  each  trade  study  and  provides  an  analysis  of 
each  trade  study.  Chapter  Five  explains  conclusions  and  recommendations  for  further 
research. 
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2  Background 


2.1  Introduction 

This  chapter  describes  the  fundamental  topics  that  became  the  building  blocks  of 
this  research.  There  is  a  brief  overview  of  GPS  along  with  an  explanation  of  system 
architecture  and  pseudorange  measurements.  This  chapter  describes  the  applicable 
reference  frame  coordinate  systems,  explains  the  theory  and  mechanics  of  Two-Way 
Time  Transfer,  and  demonstrates  the  method  of  using  a  Kalman  filter. 


2.2  GPS  Overview 

The  following  sections  will  describe  the  GPS  architecture,  including  the 
Operational  Control  Segment  (OCS),  the  space  segment,  and  the  user  segment  as  shown 
in  Figure  2.1.  The  next  section  will  describe  how  the  system  produces  navigation 
information  via  pseudorange  measurements.  This  section  closely  follows  the  descriptions 
given  in  [8]  and  [25]. 

2. 2. 1  Operational  Control  Segment  ( OCS) 

The  OCS  is  comprised  of  the  Master  Control  Station  (MCS),  monitor  stations, 
and  ground  antennas.  The  MCS  provides  the  command  and  control  for  the  entire  GPS 
system  and  tracks  GPS  satellite  orbits,  monitors  and  sustains  GPS  satellite  health,  and 
maintains  GPS  time  [16].  Operated  by  the  MCS,  the  monitor  stations  are  spread  across 
the  Earth  and  constantly  observe  GPS  satellite  signals.  Monitor  stations  have  a  GPS 
receiver,  an  atomic  clock,  meteorological  equipment,  and  communications  hardware  for 
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transmitting  measurement  data  back  to  the  MCS  [16].  The  ground  antennas  reside  with  a 
monitor  station  and  communicate  with  the  GPS  satellites.  Operated  by  the  MCS,  the 
ground  antennas  are  responsible  for  transmitting  commands  and  other  messages  to  the 
satellites  and  for  receiving  telemetry  and  system  health  from  the  satellites  [16]. 


THE  MAJOR  SEGMENTS  OF  THE  GPS  SYSTEM 


RF  SIGNAL  LEVELS 

(MMMUM  VALUES  RECEIVED) 


PLANNED  CONSTELLATION:  ^ - 

6  PLANES 

3  satellites  each  plane/ space  SEGMENT 
55®  INCLINATION  / 


UPLINK  DATA 

•  SATELLITE  EPHEMERIS 
(POSITION)  CONSTANTS 

•  CLOCK-CORRECTION 
FACTORS 

•  ATMOSPHERIC  DATA 

•  ALMANAC 


DOWNLINK  DATA 
•  CODED  RANGING  SIGNALS 
.  POSITION  INFORMATION 
.  ATMOSPHERIC  DATA 
.  ALMANAC  ^ - 


MONITOR 

STATIONS 


USER  SEGMENT 


MASTER  CONTROL 
STATION  / 


GROUNDV 

ANTENNA' 


CONTROL  SEGMENT 


TRANSMISSION  BANDS 

DATA  LINK 

FREQUENCY 

(MKi) 

L,  (OOWMUNK) 

Lj  (DOWNLINK) 

1227.6 

S-MNO  (OOWMUMQ 

2227.5 

S-AAND  (UAUNK) 

Figure  2.1:  Major  Segments  of  the  GPS  System  [15] 


2. 2. 2  Space  Segment 

The  nominal  GPS  constellation  consists  of  24  satellites  placed  into  six  orbital 
planes  at  an  altitude  of  20,200  kilometers.  Each  orbital  plane  has  an  inclination  of  55 
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degrees  and  holds  four  satellites.  The  orbital  planes  are  labeled  A  through  F,  as  seen  in 
Figure  2.2. 


D 


Figure  2.2:  The  GPS  Satellite  Constellation  [12] 


The  responsibilities  of  the  OCS  include  creating  an  ephemeris  for  each  satellite, 
which  is  a  compiled  set  of  state  vectors  for  each  satellite  predicted  over  time  [23]. 
Ephemeris  values  describe  the  shape  of  the  satellite’s  orbit,  and  can  be  used  to  calculate 
the  satellite’s  position  and  velocity.  A  Kalman  filter  calculates  ephemeris  values  by 
propagating  the  GPS  satellites’  positions  and  velocities  forward  in  time.  The  final  result 
is  called  the  broadcast  ephemeris.  All  satellites  in  the  GPS  constellation  transmit  the 
broadcast  ephemeris  data  for  each  satellite,  allowing  a  user  to  know  the  position  of  each 
GPS  satellite. 

A  precise  ephemeris  data  set  is  compiled  roughly  two  weeks  after  the  actual  orbit 
by  a  private  conglomeration  of  users.  Technically  it  is  not  a  true  ephemeris,  but  rather  a 
historical  collection  of  observables.  It  includes  measurement  data  from  hundreds  of 
ground  stations  and  accounts  for  ionospheric,  tropospheric,  and  many  other  errors.  The 
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broadcast  ephemeris  is  a  prediction,  while  the  precise  ephemeris  is  compiled  using  past 
measurements.  The  precise  ephemeris  is  extremely  accurate,  though  post-dated,  and  is 
considered  to  be  the  true  position  of  the  GPS  satellite.  Precise  ephemeris  data  used  in 
this  simulation  comes  from  the  International  GPS  Service  (IGS). 

Each  GPS  satellite  broadcasts  a  global  navigational  message  on  two  L-band 
frequencies,  LI  (1575.42  Hz)  and  L2  (1227.60  Hz)  [18].  The  navigational  message  is 
unique  to  each  satellite,  modulated  on  a  pseudo-random  noise  (PRN)  signal.  It  provides 
satellite  positions  and  satellite  clock  corrections  for  use  in  computing  the  user’s 
navigation  solution,  to  include  position,  velocity,  and  time. 

2.2.3  User  Segment 

The  GPS  signal  is  globally  available  for  use  by  any  receiver.  A  typical  GPS 
receiver  must  have  seven  components  in  order  to  correctly  acquire  and  track  the  GPS 
signal  [18],  It  must  have: 

•  an  omnidirectional  antenna  to  receive  the  encoded  signals  transmitted  by 
the  GPS  satellites 

•  a  filter  to  remove  interfering  signals 

•  an  amplifier  to  boost  the  GPS  signal 

•  a  delay  lock  loop  receiver  and  demodulator  to  provide  estimates  of  the 
pseudorange,  carrier  phase,  and  navigational  data  for  each  satellite 

•  a  navigation  data  processor  to  calculate  the  position  of  each  satellite  based 
on  the  navigation  data. 

•  an  algorithm  to  estimate  the  user  position  and  velocity  state  vector 
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•  a  reference  oscillator  to  provide  a  time  and  frequency  reference  for  the 
receiver. 

2.2.4  Pseudorange  Measurements 

Two  pseudorandom  noise  (PRN)  codes,  the  Coarse-Acquisition  (C/A)  code  and 
the  Precision  (P(Y))  code,  are  modulated  onto  the  LI  and  L1/L2  bands  respectively  [23]. 
Each  PRN  code  is  unique  for  a  given  satellite.  To  make  a  GPS  measurement,  the  PRN 
signal  transmit  time  between  the  GPS  satellite  and  the  user  must  be  determined.  To  do 
this,  the  user’s  receiver  compares  an  internal  copy  of  the  PRN  code  with  the  PRN  signal 
received  from  the  satellite.  By  calculating  the  time  shift  required  to  align  the  two  PRN 
codes  and  multiplying  by  the  speed  of  light,  the  user  can  determine  the  pseudorange 
between  the  receiver  and  the  GPS  satellite. 

Because  the  true  range  from  the  receiver  to  the  satellite  is  corrupted  by  receiver 
clock  errors,  it  is  called  a  pseudomngQ.  A  pseudorange  measurement  (p)  can  be 
expressed  as: 

p = pp'  -  f + O'"  -  pj  +  f  +  cStm  -  esc + U„  (2.1) 

where 

xs<" ,  ysat ,  zsat  =  true  ECEF  position  of  the  satellite  (meters) 
xrec,yrec,zrec  =  true  ECEF  position  of  the  receiver  (meters) 

5trec  =  receiver  clock  bias  (seconds) 

8tsat  =  satellite  clock  bias  (seconds) 

vPR  =  pseudorange  error  (meters) 

c=  speed  of  light  (2.99792458  x  108  meters/second) 
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Figure  2.3  displays  the  concept  of  a  pseudorange  measurement.  At  least  four 
satellites  are  needed  to  simultaneously  estimate  the  receiver’s  position  and  clock  error  (x, 
y,  z,  dtrec).  In  Figure  2.3,  the  b  term  is  equal  to  -cStrec  +  cStsat  - oPR ,  as  described  above. 


()tk\  y^k\  zM) 


p<  i) 

{p^}:  Pseudoranges  (measurements) 

{(xM,  yM,  z(fc))}:  Satellite  positions  (known) 

p<*>  =  ^(xw-x)2  +  (y<*>-y)2  +  (zl^-z)2  -  b 
k=  1,2 . K 

It  K  >  4,  solve  for  user  position  (x,  y,  z), 
and  receiver  clock  bias  b 


Figure  2.3:  Illustration  of  Pseudorange  Measurement  [16] 


2.3  Reference  Coordinate  Systems 

All  motion  in  the  universe  is  relative,  so  there  are  many  different  ways  to  express 
the  reference  frame  that  describes  an  object’s  orbital  motion  around  the  Earth.  Each 
reference  frame  has  advantages  for  particular  circumstances,  so  this  research  utilizes  two 
reference  frames:  Earth-Centered  Inertial  (ECI)  and  Earth-Centered  Earth-Fixed  (ECEF). 
Both  the  ECI  and  ECEF  reference  frames  are  Cartesian  (orthogonal)  coordinate  systems. 
This  section  closely  follows  the  descriptions  given  in  [25]. 
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2.3.1  ECI  Reference  Frame 


The  Earth-Centered  Inertial  coordinate  system  is  typically  used  when  describing 
the  motion  of  an  object  orbiting  the  Earth,  as  it  is  an  inertial  frame  and  is  not  concerned 
with  the  rotation  of  the  Earth.  For  ECI,  the  origin  is  located  at  the  center  of  the  Earth 

with  the  X  axis  pointing  toward  the  vernal  equinox,  also  called  the  First  Point  of  Aries,  in 
the  equatorial  plane.  TheZ  axis  points  through  the  North  Pole  along  the  Earth’s  axis  of 
rotation,  and  the  Y  axis  points  ninety  degrees  from  the  X  axis  in  the  equatorial  plane, 
following  the  right-hand  rule.  The X, Y, Z axes  are  often  referred  to  as/,  J, K  ,  as  seen  in 

Figure  2.4.  Since  the  A  axis  direction  is  fixed  in  space,  the  ECI  coordinate  system  does 
not  change  as  the  Earth  spins  about  its  axis  and  revolves  around  the  sun. 
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2.3.2  ECEF  Reference  Frame 


The  Earth-Centered  Earth-Fixed  coordinate  system  is  typically  used  for  Earth- 
based  satellite  tracking  operations  where  the  satellite’s  coordinates  need  to  be  expressed 
in  relation  to  a  geographic  location  on  the  surface  of  the  Earth.  Unlike  the  ECI  reference 
frame,  the  ECEF  reference  frame  is  always  aligned  with  a  particular  longitude,  normally 
the  Prime  Meridian  at  Greenwich,  and  thus  rotates  with  the  Earth.  The  origin  is  still 

located  at  the  center  of  the  Earth,  and  theZ  axis  still  points  through  the  North  Pole  along 
the  axis  of  rotation,  while  the  X  axis  points  toward  a  designated  longitude  and  rotates 

with  the  Earth.  The  Y  axis  always  points  ninety  degrees  from  the  X  axis.  Figure  2.5 
portrays  the  ECEF  coordinate  system. 


The  ECEF  coordinate  system  proves  useful  when  examining  geostationary  orbits, 
since  a  geostationary  orbit  rotates  at  the  same  rate  as  the  Earth.  If  a  satellite  were  situated 

in  a  perfectly  geostationary  orbit,  the  satellite’s X, Y, Z coordinates  would  remain 
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constant.  Because  this  research  examines  the  errors  in  the  navigation  solution  of  a 
geostationary  satellite,  the  results  will  be  presented  in  the  ECEF  coordinate  system, 
making  them  more  intuitive  to  interpret. 


2.4  Two-Way  Time  Transfer  Overview 

Two-Way  Time  Transfer  (TWTT)  is  a  technique  in  which  signals  are 
simultaneously  exchanged  between  two  users  to  measure  their  relative  clock  offsets.  If 
the  paths  between  the  two  users  are  reciprocal,  the  delays  cancel  and  the  difference 
between  the  two  clocks  is  half  the  difference  in  time  interval  counter  readings  [13]. 
Figure  2.6  explains  the  TWTT  technique  in  more  detail  using  a  simplified  setup  with  two 
clocks  connected  via  two  cables. 


Clock  A  Clock  B 


Figure  2.6:  Two-Way  Time  Transfer  Technique  [13] 


In  Figure  2.6, 

A  =  Time  Interval  Counter  measurement  from  Clock  A 
B  =  Time  Interval  Counter  measurement  from  Clock  B 
dAB  =  delay  caused  in  cable  from  Clock  A  to  Clock  B 
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cIba  =  delay  caused  in  cable  from  Clock  B  to  Clock  A 

R(A)  =  difference  between  Clock  A  and  Clock  B  with  cable  delay  deA 

R(B)  =  difference  between  Clock  B  and  Clock  A  with  cable  delay  dAB 

There  are  two  types  of  TWTT,  those  being  static  and  dynamic.  Static  TWTT  uses 
two  or  more  transceivers  whose  positions  are  held  constant  during  the  transmission  and 
reception  of  the  measurement  signals.  Dynamic  TWTT  is  a  more  recent  development 
that  allows  one  or  more  of  the  transceivers  to  be  moving  [2],  This  section  closely  follows 
the  descriptions  given  in  [8]  and  [25]. 

2.4.1  Static  TWTT 

TWTT  commonly  involves  two  static  clocks  on  the  Earth’s  surface  that  use  a 
geostationary  satellite  to  relay  signals  between  them.  The  connection  between  the  two 
clocks  proceeds  through  an  antenna,  a  transmitter,  an  uplink  to  the  satellite,  a  route 
through  the  satellite,  a  downlink  from  the  satellite,  an  antenna,  and  a  receiver  [13].  This 
arrangement  is  depicted  in  Figure  2.7. 


Figure  2.7:  Static  Two-Way  Time  Transfer  Using  a  Satellite  [13] 
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In  Figure  2.7, 

dAs  =  delay  between  Receiver  A  and  satellite  during  time  of  transmission 
dsA  =  delay  between  satellite  and  Receiver  A  during  time  of  transmission 
cIbs  =  delay  between  Receiver  B  and  satellite  during  time  of  transmission 
dsB  =  delay  between  satellite  and  Receiver  B  during  time  of  transmission 
diA  =  delay  in  Transmitter  A 
diB  =  delay  in  Transmitter  B 
dRA  =  delay  in  Receiver  A 
dRB  =  delay  in  Receiver  B 

dsAB  =  delay  in  satellite  when  signal  is  going  from  Transmitter  A  to 
Receiver  B 

dsBA  =  delay  in  satellite  when  signal  is  going  from  Transmitter  B  to 
Receiver  A 

TIC  =  Time  Interval  Counter 


The  Time  Interval  Counters  (TICs)  make  basic  time  interval  measurements  at 
each  site.  Each  TIC  starts  with  a  pulse  from  the  local  clock  and  ends  when  the  signal  is 
received  from  the  remote  clock.  The  remote  clock  sends  a  pulse  at  the  same  time  that  the 
local  clock  sends  a  pulse.  Nominally,  one  pulse  per  second  is  the  rate  used  for  TWTT. 
Each  station  records  the  time  interval  and  the  information  is  sent  to  the  other  station 
where  the  two  values  are  differenced.  The  data  rate  required  for  exchanging  pulses  is 
trivial,  allowing  TWTT  to  occur  in  real  time  [13]. 

The  time  interval  information  that  is  differenced  at  each  station  includes  all 
delays  shown  in  Figure  2.7.  An  additional  delay  term,  the  Sagnac  delay,  is  also  included. 
The  Sagnac  delay  is  associated  with  the  Earth’s  rotation  and  the  fact  that  transmitted 
signals  have  a  finite  velocity,  as  demonstrated  in  Figure  2.8.  The  Earth  stations  and 
satellite  are  at  position  1  when  the  pulses  are  sent  to  the  satellite.  The  rotation  of  the 
Earth  causes  the  stations  to  move  to  position  3  before  the  signals  are  received.  The 
Earth’s  rotation  and  the  finite  velocity  of  the  signals  have  combined  to  make  the  transmit 
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signal  path  longer  from  Station  B  to  Station  A,  and  shorter  from  Station  A  to  Station  B. 
For  static  TWTT,  the  Sagnac  error  is  effectively  a  deterministic  propagation  delay. 


Satellite  Motion 


Figure  2.8:  Demonstration  of  the  Sagnac  Delay:  Earth’s  Rotation  Produces  Non-Reciprocity  [13] 
By  adding  all  delays  to  the  difference,  the  complete  time  interval  measurement 

for  each  station  can  be  stated  in  the  following  equations: 


R(A)  -  A  -  B  +  diB  +  des  +  dsBA  +  dsA  +  dRA  +  Sab  (2.2) 

R(B)  =  B  -  A  +  diA  +  dAs  +  dsAB  +  dsB  +  dRB  +  Sba  (2.3) 


where 

R(A)  =  time  interval  counter  reading  for  Station  A 
R(B)  =  time  interval  counter  reading  for  Station  B 
Sab  =  Sagnac  delay  from  Station  A  to  Station  B 
Sba  =  Sagnac  delay  from  Station  B  to  Station  A 
A  =  time  of  Clock  A 
B  =  time  of  Clock  B 

All  other  variables  are  the  same  as  Figure  2.7 
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Differencing  Equations  2.2  and  2.3,  as  seen  in  Figure  2.9,  produces  the  expanded 


Two-Way  Time  Transfer  equation. 


A-B-  -RW  '  R<B> 

2 

T\ me  Interval 
Counter  Readings 

dTB  dRB  dTA  dRA 

2  +  2 

Earth  Station 
Equipment 

dASdSA  dBS‘dSB 

2  2 

Propagation 

dSABdSBA 

2 

Delay  in  Satellite 

2wAr 

C2 

Earth's  Rotation 

Figure  2.9:  The  Two-Way  Time  Transfer  Equation  [13] 

In  Figure  2.9, 

C  =  speed  of  light 

Ar  =  area  enclosing  the  projection  of  the  satellite  onto  the  Earth’s 
equatorial  plane 
co  =  Earth  rotation  rate 


In  the  case  of  static  TWTT,  the  uplink  and  downlink  propagation  delays  are 
essentially  equal,  thus  dsA  ~  dAs  and  dsB  ~  dss-  The  difference  in  Sagnac  delays  (Sab  - 
Sba)  is  constant  for  the  static  case  of  TWTT,  and  the  transmitting  and  receiving  delay  for 
the  Earth  station  equipment  is  equal,  as  is  the  satellite  delay  when  routing  a  signal  from 
Station  A  to  Station  B  or  Station  B  to  Station  A.  These  factors  cause  cancellation  when 
the  measurements  are  differenced,  and  the  TWTT  equation  is  reduced  to: 

A  -  B  =  i  [R(A)  -  R(B)  +  Sab  -  Sba]  (2.4) 

where 

R(A)  =  time  interval  counter  reading  for  Station  A 
R(B)  =  time  interval  counter  reading  for  Station  B 
Sab  =  Sagnac  delay  from  Station  A  to  Station  B 
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Sba  =  Sagnac  delay  from  Station  B  to  Station  A 
A  =  time  of  Clock  A 
B  =  time  of  Clock  B 

Great  precision  and  accuracy  can  be  achieved  by  using  a  geostationary  satellite 
for  comparing  and  synchronizing  clocks  in  this  manner.  This  accuracy  is  a  result  of 
many  simplifications  that  occur  in  the  process  of  calculating  a  time  difference 
measurement  using  TWTT.  Propagation  delays  that  occur  during  uplink  and  downlink  to 
and  from  the  satellite  are  essentially  equivocal,  causing  them  to  cancel.  Sagnac  effects 
can  be  accurately  calculated  without  the  need  for  precise  information  on  clock  locations 
and  errors,  further  reducing  complexity. 

2.4.2  Dynamic  TWTT 

Dynamic  TWTT  is  accomplished  in  the  same  fashion  as  static  TWTT  with  the 
exception  that  one  or  more  of  the  receivers  is  moving.  The  moving  receiver(s)  introduce 
motion-related  errors  that  were  not  present  in  the  static  case,  changing  the  TWTT  time 
differencing  equation.  A  dynamic  TWTT  configuration  is  illustrated  in  Figure  2.10. 


Figure  2.10:  Dynamic  Two-Way  Time  Transfer  Using  a  Satellite  [3] 
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Not  all  of  the  cancellations  that  applied  to  the  static  TWTT  case  transfer  to  the 
dynamic  case.  For  the  dynamic  TWTT  scenario  in  Figure  2.10,  it  can  be  assumed  that 
dAs  ~  dsA,  since  the  geostationary  satellite  has  no  relative  motion  with  respect  to  the  earth 
station  and  the  path  length  does  not  change.  This  replicates  the  situation  in  the  static 
case.  Unlike  the  static  case,  dsB  ^  des  for  the  dynamic  case,  since  the  mobile  platform 
has  moved  during  the  transmission  of  signals,  causing  the  transmit  and  receive  path 
lengths  to  be  different  between  the  geostationary  satellite  and  mobile  platform.  Because 
the  mobile  platform  is  in  motion,  the  Sagnac  effect  will  also  vary  and  produces  a  time- 
dependent  value. 

Taking  all  this  into  account,  the  time  differenced  measurement  for  dynamic 
TWTT  becomes: 

A  -  B  =  ^  [R(A)  -  R(B)  -  A Propagation  Delay  +  Sab  -  Sba]  (2.5) 

where 

R(A)  =  time  interval  counter  reading  for  Station  A 
R(B)  =  time  interval  counter  reading  for  Station  B 
Sab  =  Sagnac  delay  from  Station  A  to  Station  B 
Sba  =  Sagnac  delay  from  Station  B  to  Station  A 
A  =  time  of  Clock  A 
B  =  time  of  Clock  B 

A Propagation  Delay  =  change  in  propagation  delay  over  measurement 
interval 

The  A Propagation  Delay  term  is  a  time-varying  value  that  changes  based  on  the 
relative  motion  of  the  mobile  platform  as  well  as  how  the  velocity  vector  is  projected 
onto  the  line  of  sight  vector  from  the  geostationary  satellite.  The  Sagnac  delay  term  (Sab 
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-  Sba)  is  also  time-varying,  changing  based  on  the  absolute  position  of  the  two  receivers 
and  the  velocity  vector  projected  onto  the  equatorial  plane  [3]. 

2.5  Kalman  Filter 

A  Kalman  filter  is  a  sequential  filter  that  continuously  improves  the  estimate  of  a 
state  vector  by  sequentially  incorporating  new  data  measurements  into  the  estimate  as 
they  become  available  [24],  A  state  vector  can  include  data  regarding  position,  velocity, 
acceleration,  and  any  other  pertinent  information.  A  Kalman  filter  predicts  the  state  in  a 
future  time  epoch  and  then,  upon  reaching  that  new  epoch,  it  uses  collected 
measurements  to  correct  the  estimated  state  before  predicting  again  to  a  new  epoch. 
Kalman  filters  are  able  to  compensate  for  ill-known  or  incompletely  modeled  dynamical 
systems  [23]  and  are  ideal  for  handling  stochastic  systems,  such  as  navigation.  This 
section  uses  the  equations  and  descriptions  provided  in  [22]  for  presenting  a  Kalman 
filter  summary. 

The  Kalman  filter  is  an  iteration  of  time  propagation  and  measurement 
incorporation  [22],  It  begins  with  an  initial  estimated  state  vector  x(/0)  and  an  initial 

estimated  covariance  P (t0) .  The  covariance  matrix  is  a  statistical  measure  of  the 
accuracy  of  the  state  vector  x(70).  In  the  next  step  the  state  x(/0)  and  covariance  P(t0) 

are  predicted  for  the  next  epoch.  These  predictions  are  merely  calculated  and  are  not 
computed  using  measurement  data.  The  Kalman  filter  uses  Equations  2.6  and  2.7  to 
propagate  the  state  and  the  covariance  forward  in  time  from  t+_ x  to  fk .  Using  this 
notation,  tk  and  tk[  refer  to  the  times  when  measurements  are  incorporated  and  a 
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superscript  plus  sign  (+)  indicates  the  specified  time  is  after  measurement  incorporation, 
while  a  superscript  minus  sign  (-)  indicates  the  specified  time  is  before  measurement 
incorporation. 


x(4)  =  ®(4-4-i)i(Ci) 

P  (4 )  =  ®(4  -4-.)P(4+-i)^(4  -4-i)+Q, 


(2.6) 

(2.7) 


where 

x(4_, )  =  current  state  vector  (before  propagation) 
x(4)  =  state  vector  propagated  to  next  epoch 
P(4+-  a)  =  current  covariance  matrix  (before  propagation) 
p(4"  )  =  covariance  matrix  propagated  to  next  epoch 
0(4  -  4  j )  =  state  transition  matrix,  which  propagates  the  state  and 
covariance  forward  in  time 

Qd  =  second  moment  of  the  process  noise,  which  is  the  discrete  time 

equivalent  covariance  associated  with  process  noise  that  occurs 
when  covariance  P  is  propagated  through  time  [23] 

Once  the  state  vector  and  covariance  matrix  are  propagated  forward  in  time  to  the 

next  epoch,  real  measurements  are  incorporated  to  correct  any  error  that  may  have  been 

introduced  during  propagation.  The  measurements  must  first  be  weighted  properly 

before  they  can  be  incorporated  into  the  estimated  state  and  covariance.  If  the 

measurements  are  very  poor,  the  filter  must  know  to  rely  less  on  the  collected  data  and 

depend  more  upon  its  own  estimation  of  the  state.  Alternatively,  if  the  measurements  are 

very  good,  the  filter  must  know  to  rely  more  heavily  upon  the  collected  data  and  depend 

less  upon  its  own  estimation  of  the  state.  This  weighting  factor  is  calculated  using  the 

Kalman  gain,  as  seen  in  Equation  2.8. 

K  =  P(4  )Hr[HP(4  )Hr  +  R]-1  (2.8) 


where 

K  =  Kalman  gain  matrix 
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p  it;  )  =  covariance  matrix  propagated  to  next  epoch  (but  before 
measurement  incorporation) 

H  =  observational  partial  derivative  matrix,  which  relates  the  linearized 
observations  (z)  to  the  estimated  states  (x) 

R  =  measurement  noise  covariance  matrix 


Once  the  Kalman  gain  has  been  calculated,  the  state  and  covariance  are  updated 
by  incorporating  the  measurement  data.  This  process  is  illustrated  in  Equations  2.9  and 
2.10. 


x(4+  )  =  x(t~  )  +  K[z  -  Hx(4  )]  (2.9) 

P(4+)  =  P(4)-KHP(4)  (2.10) 


where 

z  =  measurement  vector 

z-Hx(t“)  =  r  =  residual  vector,  which  indicates  how  much  error  exists 

between  the  estimated  state  and  the  measurements  (want  residuals 
to  be  small,  as  it  indicates  the  estimate  is  accurate) 


At  this  point,  the  cycle  repeats  as  the  Kalman  filter  again  propagates  the  state  and 
covariance  to  the  next  epoch,  incorporates  collected  measurements,  and  updates  the  state 
and  covariance  using  the  measurements.  The  Kalman  filter  loop  is  displayed  in  Figure 
2. 1 1  (using  slightly  different  notation). 
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Enter  prior  estimate  x0  and 
its  error  covariance  P0 


t 


Figure  2.11:  Kalman  Filter  Loop  [7] 


2.6  Summary 

This  chapter  presented  the  relevant  background  information  necessary  for  this 
research.  Several  key  issues  were  explained,  including  the  GPS  system  and  pseudorange 
measurements,  Earth-based  reference  coordinate  systems,  both  static  and  dynamic  Two- 
Way  Time  Transfer  techniques,  and  the  Kalman  fdter.  Chapter  3  will  describe  the 
methodology  of  this  research  and  how  it  used  the  background  concepts  explained  in 
Chapter  2. 
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3  Methodology 


3.1  Introduction 

This  chapter  describes  the  methodologies  and  assumptions  that  were  used  to 
accomplish  the  research  objectives  described  in  Chapter  1.  This  research  is  based  on  a 
simulation  created  using  MATLAB®  and  contains  five  main  functions. 

The  first  main  function  involves  collecting  desired  input  parameters  from  the 
user.  The  second  function  uses  the  input  parameters  to  create  truth  data  that  will  simulate 
the  environment  that  is  being  measured.  The  third  function  uses  the  truth  data  to  generate 
pseudorange  and  TWTT  measurements  for  a  geostationary  satellite.  The  fourth  function 
inputs  the  generated  measurements  into  a  Kalman  filter  and  predicts  the  state  of  the 
satellite  at  each  epoch  in  the  simulation.  The  fifth  and  final  function  takes  the  results  of 
the  Kalman  filter  and  compares  them  to  the  truth  data  to  determine  the  accuracy  of  the 
filter. 


3.2  Parameters 

The  simulation  begins  by  collecting  all  the  desired  input  values  for  a  host  of 
variables  that  will  be  used  throughout  the  simulation.  The  list  of  input  variables  is 
described  in  the  following  sections.  All  the  values  in  the  parameters  function  are 
declared  globally  for  use  in  all  functions  throughout  the  simulation. 
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3.2.1  Initial  ECI  State  for  Geostationary  Satellite 


The  initial  ECI  state  consists  of  X,  Y,  Z  position  and  velocity  values  that  describe 
the  location  of  the  satellite  at  the  start  of  the  simulation.  This  vector  is  nominally 
configured  to  place  the  satellite  over  the  Prime  Meridian  (0  degrees  longitude)  at  an 
altitude  of  35,785.863  kilometers  with  an  orbital  velocity  of  3.07466  km/s.  The  result  is 
a  circular  geostationary  orbit. 

3.2.2  Simulation  Run  Time  and  Time  Step  Interval 

The  simulation  run  time  is  set  for  a  single  day.  The  time  step  interval  denotes 
how  often  measurement  data  is  collected  by  the  satellite  as  it  orbits  and  is  set  for  60 
seconds.  If  a  higher  measurement  resolution  is  desired,  the  time  step  can  be  shortened  at 
the  expense  of  a  longer  computation  time. 

3 .2.3  Ephemeris  Date  Selection 

As  discussed  in  Chapter  2,  the  ephemeris  describes  the  location  of  a  satellite 
within  its  orbit.  The  ephemeris  data  for  the  GPS  constellation  may  be  chosen  from  any  of 
ten  possible  days  that  were  randomly  selected  from  each  year  between  1997  and  2006. 
For  simplicity  in  the  model,  each  chosen  day  is  a  Sunday.  GPS  time  is  recorded  in  week 
seconds  and  resets  at  midnight  every  Sunday.  By  starting  the  simulation  at  midnight  on 
Sunday,  which  is  zero  GPS  week  seconds,  many  complications  can  be  avoided  in  the 
simulation  code.  These  issues  arise  when  attempting  to  keep  track  of  the  GPS  week 
second  time  interval  and  initialize  the  time  vector  accordingly,  making  it  much  easier  to 
always  start  at  a  time  value  of  zero. 
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3.2.4  GPS  Measurement  Model  Selection 


Two  models  are  available  for  generating  GPS  pseudorange  measurements — a 
simple  model  and  an  advanced  model.  The  simple  model  uses  a  constant  cut-off  angle 
between  the  GPS  satellite  and  the  Earth.  As  long  as  the  Earth  is  not  blocking  the 
transmitted  signal  from  the  GPS  satellite,  it  is  received  by  the  GEO  satellite.  The 
advanced  model  is  more  realistic,  because  it  accounts  for  the  antenna  gain  on  the  GPS 
satellite  and  the  GEO  satellite  and  includes  the  signal  path  loss  that  occurs  as  the  signal 
travels  through  space.  If  the  received  power  is  below  a  specified  threshold,  the  signal  is 
considered  too  weak  and  cannot  be  used  for  generating  pseudorange  measurements. 

The  advanced  GPS  measurement  model  has  additional  parameters  that  can  be 
adjusted.  There  are  two  gain  patterns  that  can  be  used  for  the  GPS  satellite,  simple  and 
advanced.  The  simple  gain  pattern  contains  fewer  values,  ranging  between  0-30  degrees 
off-boresight,  and  is  based  on  data  found  in  [19].  The  advanced  gain  pattern  has  values 
from  0-180  degrees  off-boresight  and  more  realistically  portrays  an  actual  GPS  antenna 
and  is  based  on  data  found  in  [9].  The  GEO  satellite  uses  the  gain  pattern  taken  from  the 
commercial  patch  antenna  listed  in  [19],  replicating  actual  flight  hardware. 

3.2.5  Clock  Type  Selection 

The  simulation  has  a  choice  of  three  clock  types  to  use  as  the  GEO  satellite  clock 
and  as  the  TWTT  reference  clock.  The  possible  clock  oscillators  include  Rubidium, 
Cesium,  and  ovenized  crystal.  The  clock  parameters  are  taken  from  [14].  The  Rubidium 
clock  is  the  most  precise  and  the  ovenized  crystal  is  the  least  precise.  Ovenized  crystal  is 
an  inexpensive  and  prolific  type  of  clock  oscillator  and  is  the  typical  clock  of  choice  built 
into  current  satellites. 
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3.2.6  Clock  Model  Parameters 


In  order  to  accurately  predict  how  the  clocks  will  propagate  forward  in  time, 
certain  parameters  must  be  configured.  The  GEO  satellite  clock  and  TWTT  reference 
clock  are  assigned  an  initial  accuracy  in  bias  and  drift  that  is  fed  into  the  truth  state 
propagator  and  the  Kalman  filter,  allowing  both  propagators  to  model  the  clocks 
accurately.  These  initial  values  are  standard  deviations  of  the  bias  and  drift  error  and  will 
be  multiplied  by  a  random  number  and  fed  into  the  covariance  matrix  when  starting  the 
clock  simulation. 

3.2. 7  Two-  Way  Time  Transfer  Parameters 

The  TWTT  measurements  are  modeled  to  have  a  white  Gaussian  measurement 
noise  error,  expressed  as  a  standard  deviation.  The  TWTT  measurement  noise  error  is  a 
result  of  signal  degradation  during  transmission  and  in  the  hardware,  affecting 
measurement  accuracy,  and  will  be  held  constant.  The  values  used  were  10,  3,  0.3,  0.03, 
and  0.003  meters  for  this  simulation. 

3.2.8  Kalman  Filter  Parameters 

The  Kalman  filter  contains  initialization  values  that  are  used  in  the  covariance 
matrix  for  predicting  the  level  of  measurement  accuracy.  The  initial  accuracy  in  the 
position  and  velocity  measurements  are  provided  to  the  filter  when  it  starts.  There  is  also 
a  process  noise  value,  expressed  as  a  standard  deviation,  which  is  a  degradation  of  the 
measurement  signal  as  it  travels  through  the  associated  hardware,  reducing  accuracy. 
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3.2.9  On/Off  and  Selection  Flags 


Numerous  naturally-occurring  errors  are  built  into  the  simulation  model  that  can 
be  turned  on  or  off.  The  on/off  flags  enable  or  disable  clock  errors  in  the  GEO  satellite 
and/or  TWTT  reference  clock,  measurement  noise  errors  in  the  pseudorange  and/or 
TWTT  measurements,  hardware-related  process  noise  errors,  and  initial  state  errors  in  the 
Kalman  filter. 

The  selection  flags  allow  the  user  to  choose  different  modeling  options  in  the 
simulation.  The  GPS  measurement  model  includes  a  simple  and  advanced  model,  and 
the  TWTT  measurement  model  includes  two  choices:  using  no  TWTT  measurements 
and  using  TWTT  measurements  where  the  reference  clock  has  a  drift  rate  that  does  not 
coincide  with  GPS  time. 

3.2.10  Monte  Carlo  Parameters 

The  simulation  is  equipped  to  handle  Monte  Carlo  simulation,  involving  multiple 
iterations  of  the  simulation  with  different  random  numbers.  Monte  Carlo  simulation  is  a 
statistical  analysis  that  strives  to  provide  enough  data,  through  repetition  and  random 
inputs,  to  accurately  represent  the  stochastic  characteristics  of  the  model.  In  this 
simulation,  Monte  Carlo  simulation  is  turned  on  or  off  with  a  flag,  and  the  total  number 
of  iterations  is  defined  by  the  user. 

3.2.11  Constants 

This  simulation  uses  many  constant  values  throughout  the  various  sub-functions. 
These  constants  include  the  speed  of  light,  the  radius  of  the  Earth,  the  radius  of  a  GPS 
orbit,  the  frequency  of  the  GPS  LI  signal,  nominal  transmit  power  of  a  GPS  antenna,  the 
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gravitation  parameter  of  the  Earth,  and  the  rotation  rate  of  the  Earth.  To  eliminate  chance 
for  error,  the  sub-functions  reference  all  constants  from  a  single  global  location,  as 
opposed  to  being  declared  in  each  sub-function  that  requires  the  value. 


3.3  Truth  Model 

The  truth  model  function  is  responsible  for  generating  all  data  that  will  be 
considered  as  the  absolute  truth.  There  are  three  main  sections  in  this  function: 
propagating  the  geostationary  satellite  state  forward  over  a  specified  time  interval, 
calculating  the  positions  and  clock  states  for  each  GPS  satellite  during  that  time  interval, 
and  modeling  clock  errors. 


3.3.1  Propagate  Geostationary  Satellite  State 

Implementation  of  a  simple  Kalman  filter  propagates  the  GEO  satellite  state 
vector  into  the  future.  The  initial  state  vector  x(/0)  is  provided  along  with  an  initial 

covariance  matrix  P(t0)  and  a  dynamics  matrix  F(r0  ) .  As  explained  previously,  the 

covariance  matrix  describes  the  accuracy  of  the  state  vector  values,  and  the  dynamics 
matrix  explains  the  motion  of  the  state  vector. 

The  GEO  satellite  state  vector  contains  six  values,  displayed  in  Equation  3.1, 
where  X,Y,Z  are  the  position  values  and  X,Y,Z  are  the  velocity  values. 


X 

Y 
Z 
X 

Y 
Z 


(3.1) 
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The  dynamics  matrix  describes  the  motion  of  the  GEO  satellite,  which  is  a 
circular  orbit,  and  is  in  a  continuous  form.  The  dynamics  matrix  is  constructed  to  satisfy 
the  relationship  described  in  Equation  3.2. 


*  =  F*  (3-2) 

The  orbital  acceleration  for  a  two-body  circular  orbit  is  calculated  by  using 
Equation  3.3.  This  expression  is  needed  to  find  the  time  derivative  of  the  orbital  velocity, 
which  is  needed  to  build  the  relationship  indicated  in  Equation  3.2. 


M  ~ 
a  =  — x 

r 


(3.3) 


where 


a  =  orbital  acceleration  vector 
/i  =  Earth  gravitational  constant 

r  =  +72  +Z2  =  orbital  radius 

x  =  orbital  position  vector 

Using  Equation  3.3,  a  list  of  relationships  is  formulated  for  use  in  creating  the 
dynamics  matrix.  The  list  of  relationships  is  shown  in  Equation  3.4,  where 

[*  *  *  ”1  T  Y 

X  Y  Z  X  Y  Z  =  [xj  x2  x3  x4  x5  x6  ]  . 


X;  =X4 


x2=xs 


x}=x6 

*4=5*1 


*5  =fX2 


*6  =^*3 


(3.4) 


Using  Equation  3.4,  Equation  3.2  can  be  written  in  expanded  form  as: 
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0  0  0  0  1  0 

0  0  0  0  0  1 

-4  0  0  0  0  0 

r 

o  -4  o  ooo 

r 

o  0-4000 

r 

For  ease  in  a  computer  simulation,  the  continuous  form  equation  must  be 
converted  to  a  discrete-time  equation  [22],  This  conversion  takes  the  dynamics  matrix  F 
and  turns  it  into  the  state  transition  matrix  <D,  as  shown  in  Equation  3.6. 

O  =  eF/U  (3.6) 

where 

A t  =  4+1  -  tk  =  time  interval 

A  similar  conversion  must  take  place  for  the  covariance  matrix  Q  to  ensure  it  is 
also  in  a  discrete-time  form.  Shown  below,  the  continuous  Q  matrix  represents  process 
noise  being  added  to  the  acceleration  terms  to  account  for  unmodeled  orbit  perturbations. 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  (Tn2  0  0 

0  0  0  0  <7n2  0 

0  0  0  0  0  (Tn2 

(3.7) 

Using  Equation  3.8,  the  continuous  Q  becomes  the  discrete-time  Q(i.  The  Qd 
matrix  will  be  used  when  introducing  process  noise. 

At 

Qd  =  J®(r)QOT(r)Jr  (3.8) 
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The  next  step  is  to  propagate  the  state  to  the  next  time  interval  by  using  Equation 
3.9. 

x(4~ )  =  ®  (4  -  4-i  )x(4+-i )  (3-9) 

If  process  noise  is  to  be  included  in  the  propagation,  correlated  randomized  white 
noise  must  be  introduced  to  Equation  3.9.  This  is  achieved  by  using  the  U-D 
factorization  algorithm  (explained  in  [7])  on  the  Qd  matrix,  which  will  provide  the 
correlation  factor.  The  decomposition  of  Qd  into  U  and  D  is  illustrated  in  Equation  3.10. 

Qd  =  UDUt  (3.10) 

where 

U  =  composed  of  elements  along  major  diagonal  of  Qd,  the  nontrivial 
elements  in  upper  triangular  part  of  Qd,  zeros  elsewhere 
D  =  composed  of  elements  in  major  diagonal  Qd,  zeros  elsewhere 

By  multiplying  the  U-D  factorization  results  with  a  6  x  1  vector  of  random 
numbers  generated  by  MATLAB®,  the  result  is  a  6  x  1  vector  Wd  of  correlated 
randomized  white  noise  values  that  correspond  with  the  positions  and  velocites  in  the 
state  vector  x .  This  process  is  described  in  Equation  3.11. 

wd=uVDn  (3.11) 

where 

wd  =  6  x  1  vector  of  correlated  white  noise 
n  =  6  x  1  vector  of  random  numbers 

If  so  desired,  process  noise  can  be  included  in  the  state  propagation  by  using 
Equation  3.12.  Inclusion  of  process  noise  is  selected  by  using  a  flag  in  the  parameters 
function.  It  was  always  selected,  except  for  initial  simulation  debugging. 

x(4 )  =  °(4  -  4-i  )x(4+-i )  +  wd  (3  • 1 2) 
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3.3.2  Calculate  GPS  Satellite  Positions  and  Clock  States 


After  propagating  the  GEO  satellite  state  vector  forward  in  time  over  the  entire 
simulation  time  interval,  the  truth  model  function  then  calculates  the  positions  and  clock 
states  of  each  individual  GPS  satellite  over  the  entire  simulation  time  interval.  This 
information  is  stored  for  use  in  the  generated  measurements  function. 

There  are  two  sets  of  ephemeris  data,  broadcast  and  precise.  As  described  in 
Chapter  2,  the  broadcast  ephemeris  is  a  prediction  of  the  orbital  parameters  for  the  GPS 
satellites,  while  the  precise  ephemeris  is  a  calculation  of  the  orbital  parameters  using 
collected  measurement  data.  The  broadcast  ephemeris  is  used  by  a  GPS  receiver  to 
calculate  a  navigation  solution,  while  the  precise  ephemeris  contains  two-week-old  data 
and  can  only  be  used  for  post-processing.  In  this  simulation,  the  precise  ephemeris  is 
used  to  calculate  the  true  position  and  clock  state  of  each  GPS  satellites,  while  the 
broadcast  ephemeris  is  used  when  generating  simulated  measurements. 

3.3.3  Clock  Model 

GPS  depends  upon  highly  accurate  atomic  clocks  that  can  provide 
synchronization  between  the  satellites  and  receiver.  Synchronization  occurs  through 
estimating  the  time  offset,  drift,  and  drift  rate  of  the  receiver  clock  relative  to  GPS  time. 
To  ensure  the  simulation  accurately  represents  real  clocks,  the  clock  model  must  replicate 
the  performance  of  atomic  clocks.  This  sub-section  closely  follows  the  descriptions 
given  in  [8]  and  [25], 

GPS  satellites  use  Rubidium  and  Cesium  atomic  clocks  to  maintain  accurate  GPS 
time.  In  addition  to  these  clock  types,  ovenized  crystal  will  also  be  modeled  for  use  on 
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the  GEO  satellite.  Ovenized  crystal  is  a  very  common  clock  found  on  satellites  that  do 
not  require  the  extremely  precise  time  measurements  used  in  navigation. 

The  performance  of  atomic  clocks  can  be  simulated  using  a  three-state 
polynomial  process  driven  by  white  noise.  The  discrete  process  model  and  its  covariance 
can  be  written  as  [6]: 
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where 

Xj  (4 )  and  Xj  (tk+] )  =  clock  bias  error  at  times  tk  and  tk+] 
x2  (4 )  and  x2  (4+1 )  =  clock  drift  error  at  times  tk  and  tk+1 
x3  (4 )  and  x3  (4+1 )  =  clock  drift  rate  error  at  times  4  and  tk+l 
t  =  4+1  -  4  =  time  interval 

wj(k),  W2(k),  and  H’s(/c)  =  independent  white  noises 

qi,q2,q$  =  continuous  process  noise  power  spectral  densities  representing 
the  bias,  drift,  and  drift  rate 

O(r)  =  state  transition  matrix  propagating  current  clock  bias,  drift,  and 
drift  rate  errors  forward  in  time  from  4  to  tk+l 
Qd  (r)  =  discrete-time  process  noise  covariance  matrix 


The  clocks  cannot  be  modeled  deterministically  because  of  their  stochastic  nature. 
Instead,  the  performance  of  the  random  walk  noise  values  ( wi ,  W2,  ws )  is  modeled  and  the 
characteristic  Allan  Variance  curves  of  the  atomic  frequency  standards  are  matched  [7], 
An  example  of  a  three-state  random  clock  process  is  illustrated  in  Figure  3.1.  A  best  fit 
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curve  drawn  through  the  plot  shows  how  the  three-state  atomic  clock  performance  is 
quadratic  in  nature.  The  statistics  of  the  random  walk  noise  values  are  determined  by  the 
values  of  the  variance  elements  (qi,  q2,  qs)  of  Qd  in  Equation  3.15  [5]. 


The  q  values  for  this  simulation  were  chosen  by  leveraging  research  conducted  in 
the  Clock  Improvement  Initiative  [14]  and  are  displayed  in  Table  3.1.  To  calculate  a 
clock’s  three-state  random  process  in  the  simulation,  initial  clock  bias,  drift,  and  drift  rate 
values  are  selected  from  Table  3.1  and  then  propagated  using  Equation  3.14.  The  Qd 
from  Equation  3.15  was  used  to  generate  properly  correlated  wi,  W2,  and  wj  terms  using 
the  UD  factorization  technique  described  in  Equations  3.10  and  3.11. 


Table  3.1:  Process  Noise  Values  for  Atomic  Clocks 


Rubidium  Clock 

Cesium  Clock 

Ovenized  Crystal  Clock 

qi  (bias) 

1.11  x  10‘22  s2/s 

4.44  x  10‘22  s2/s 

1.6  x  10'21  s2/s 

q2  (drift) 

2.22  x  10'32  s2/s3 

3.33  x  10'32  s2/s3 

1 6n  x  10'23  s2/s3 

qs  (drift  rate) 

6.66  x  10'45  s2/s5 

0  s2/s5 

0  s2/s5 
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When  using  a  three-state  model,  a  value  of  zero  for  the  qs  term  creates  a 
singularity.  To  prevent  this,  a  value  of  1  x  lO'100  s2/s5  is  used  in  the  simulation  instead  of 
zero.  This  value  is  sufficiently  small  enough  to  represent  zero  while  also  being  large 
enough  to  prevent  a  singularity  from  occurring. 

3.4  Generated  Measurements 

The  measurement  generation  function  is  responsible  for  creating  pseudorange 
measurements  by  using  the  information  supplied  by  the  truth  generation  function.  The 
pseudorange  measurements  can  be  calculated  using  a  simple  or  an  advanced  GPS 
measurement  model. 

3. 4. 1  Pseudorange  Measurements 

Pseudorange  values  are  normalized  range  measurements  with  the  addition  of 
errors  due  to  pseudorange  measurement  noise,  GPS  satellite  clock  bias,  and  receiver 
clock  bias.  The  pseudorange  equation  is: 

P  =  -  -V )2  +  O'"  - )2  +  (z“  -  )2  +  cSt„  -  cSt -  +  u„  (3,16) 

where 

xsat ,  ysat ,  zsat  =  true  ECEF  position  of  the  satellite  (meters) 
xrec  ’  Tree  ’  Zrec =  true  ECEF  position  of  the  receiver  (meters) 

Strec  =  receiver  clock  bias  (seconds) 

Stsat  =  satellite  clock  bias  (seconds) 
oPR  =  pseudorange  error  (meters) 
c  =  speed  of  light  (meters/second) 

The  simulation  allows  various  errors  to  be  turned  on  and  off,  so  the  pseudorange 
is  calculated  in  stages  to  produce  more  efficient  code.  Initially,  the  pseudorange 
measurement  is  the  exact  distance  between  the  GPS  satellite  and  the  GEO  satellite,  as  is 
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described  by  the  square  root  value  in  Equation  3.16.  If  no  errors  were  turned  on,  this 
would  be  the  final  output  of  the  pseudorange  generation. 

The  next  step  is  to  determine  which  errors  are  turned  on,  and  add  each  error  value 
incrementally.  If  all  errors  are  turned  on  in  the  simulation,  then  the  final  output  of  the 
pseudorange  generation  would  match  Equation  3.16.  The  receiver  and  GPS  satellite 
clock  errors  are  generated  in  the  truth  model  and  would  be  incorporated  into  the 
pseudorange  measurement  if  their  on/off  flags  were  turned  on  in  the  parameters  function. 

The  measurement  generation  function  has  two  models  for  simulating  the 
pseudorange  measurements.  The  simple  model  is  very  rudimentary  and  is 
computationally  less  expensive,  while  the  advanced  model  is  much  more  realistic. 

3.4.2  Simple  GPS  Measurement  Model 

The  simple  GPS  measurement  model  uses  distinct  cut-off  angles  for  determining 
satellite  visibility  and  does  not  include  the  side  lobes  of  the  GPS  transmit  antenna.  The 
minimum  cut-off  angle  is  where  the  Earth  occludes  the  GPS  signal,  and  the  maximum 
cut-off  angle  is  where  the  GPS  antenna  no  longer  transmits.  Figure  3.2  illustrates  the 
simple  measurement  model. 

The  ionosphere  extends  to  about  400  kilometers  above  the  Earth’s  surface  and 
will  bend  the  GPS  signals,  creating  errors  in  the  pseudorange  measurement.  For  this 
reason,  the  Earth’s  surface  cannot  be  used  as  the  occlusion  for  the  GPS  signals,  since  any 
signals  that  come  close  to  the  Earth’s  surface  will  be  distorted  and  erroneous. 
Ionospheric  errors  can  be  corrected  and  removed,  but  this  simulation  does  not  use  such 
methods.  For  simplicity,  the  ionosphere  boundary  is  considered  opaque,  and  any  signals 
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that  pass  through  it  are  ignored.  In  a  sense,  it  is  like  extending  the  radius  of  the  Earth  by 
400  kilometers. 

The  nominal  beamwidth  of  a  Block  II/IIA  GPS  satellite  antenna  is  approximately 
42.6°  [1],  which  equates  to  a  21.3°  off-nadir  look  angle.  This  value  becomes  the 
maximum  look  angle  for  the  GPS  satellite,  as  seen  in  Figure  3.2. 


Figure  3.2:  Simple  GPS  Measurement  Model  Diagram  (not  to  scale) 

The  simple  GPS  measurement  model  produces  a  pseudorange  measurement  for  a 
given  GPS  satellite  only  if  the  GEO  satellite  falls  within  the  visible  region  seen  in  Figure 
3.2.  Signal  power  is  not  modeled  and  it  is  assumed  that  line-of-sight  visibility  guarantees 
a  pseudorange  measurement.  The  standard  deviation  of  the  pseudorange  measurement 
noise  error  is  constant,  set  at  five  meters.  If  the  pseudorange  measurement  noise  error  is 
not  turned  on  in  the  simulation,  the  pseudorange  value  is  the  exact  distance  between  the 
GEO  and  GPS  satellites. 
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3.4.3  Advanced  GPS  Measurement  Model 


The  advanced  GPS  measurement  model  is  much  more  robust  than  the  simple 
model,  taking  into  account  the  individual  antenna  gain  patterns  of  each  satellite  to 
determine  the  received  signal  strength.  It  then  uses  the  received  signal  strength  to 
calculate  the  value  of  the  pseudorange  measurement  noise  error,  which  is  taken  from  a 
table  of  values.  The  minimum  cut-off  angle  imposed  by  the  upper  ionosphere  is  still  in 
effect.  The  advanced  GPS  measurement  model  is  illustrated  in  Figure  3.3 


To  determine  the  strength  of  the  GPS  signal  that  is  received  by  the  GEO  satellite, 
the  satellite  nadir  look  angles  are  needed.  If  the  GEO  satellite  and  GPS  satellite  positions 
are  known,  simple  vector  math  will  produce  the  angles  9  and  a  (referenced  in  Figure  3.3), 
which  are  the  GPS  satellite  look  angle  and  GEO  satellite  look  angle,  respectively.  The 
calculated  look  angles  are  then  used  with  antenna  gain  pattern  information  to  determine 
the  received  signal  strength.  By  knowing  the  look  angle,  it  is  a  simple  measure  of  tracing 
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the  plot  to  find  the  corresponding  antenna  gain.  This  is  accomplished  for  both  the  GPS 
satellite  and  the  GEO  satellite. 

An  antenna  gain  pattern  plot  can  be  created  in  a  laboratory  by  making 
measurements  on  an  antenna  and  recording  the  signal  strength  at  various  horizontal  look 
angles.  This  simulation  uses  either  a  simple  or  advanced  GPS  antenna  gain  pattern  plot, 
both  based  on  the  Block  II/IIA  GPS  satellites.  For  the  GEO  satellite  antenna,  the 
simulation  uses  an  antenna  gain  pattern  plot  taken  from  a  commercial  patch  antenna  that 
is  currently  used  on  spacecraft. 

The  simple  GPS  satellite  antenna  gain  pattern  plot  includes  the  primary  signal 
lobe  and  the  secondary  lobe,  as  seen  in  Figure  3.4,  and  drops  off  at  approximately  43° 
from  nadir. 


The  advanced  GPS  antenna  gain  pattern  plot,  shown  in  Figure  3.5,  includes  the 
primary  and  secondary  lobes,  as  well  as  additional  information  all  the  way  out  to  180° 
from  nadir.  A  gain  pattern  is  not  necessarily  symmetrical,  but,  for  simplicity  in  this 
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simulation,  it  is  assumed  that  the  gain  pattern  does  not  change  with  varying  vertical 
angles.  Similar  GPS  broadcast  antenna  gain  pattern  data  is  located  in  [4]  and  [11]. 


Figure  3.5:  Advanced  GPS  Antenna  Gain  Pattern  Plot  [9] 


The  GEO  satellite  uses  the  gain  pattern  from  a  patch  antenna  that  flew  aboard  the 
Falcon  Gold  experiment,  designed  by  students  at  the  Air  Force  Academy  [19].  This 
particular  antenna  is  representative  of  hardware  that  has  flown  on  previous  satellites.  The 
antenna  gain  pattern  plot  is  illustrated  in  Figure  3.6. 
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Figure  3.6:  Patch  Antenna  Gain  Pattern  Plot  [19] 

The  received  power  can  be  calculated  by  using  the  transmit  and  receive  antenna 
gains,  the  transmit  power,  the  distance  between  the  transmit  and  receive  antennas,  and  the 
wavelength  of  the  transmitted  signal.  The  resulting  equation  is  listed  below  in  standard 
and  decibel  format  [16]. 


1  s,dB 


p  _  PjGjGji 

4nR2  4 n 

=  P T,dB  +  G T,dB  +  GR  dB  ~  20  log10  R  —  22  +  2AdB 


(3.17) 

(3.18) 


where 

Ps  =  received  signal  power  (Watts) 

Pt  =  signal  power  at  transmit  antenna 
Gt  =  transmit  antenna  gain 
Gr  =  receive  antenna  gain 

R  =  distance  between  transmit  and  receive  antennas 
X  =  signal  wavelength  (GPS  LI  wavelength  «  5.255  meters) 

The  resulting  value  is  the  signal  power  at  the  exit  of  the  receiver  antenna.  All 

hardware  has  an  associated  noise  value  which  must  be  subtracted  from  the  received 

power  when  the  signal  is  processed.  The  final  value  is  referred  to  as  C/No  (also  written  as 

P/No),  which  is  the  carrier  to  noise  density,  where  No  is  the  noise  density.  This 
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simulation  uses  a  standard  No  value  of  202  decibel-Watts  [16],  which  is  subtracted  from 
Equation  3.18  to  produce  the  final  C/No  value. 

Once  the  received  signal  strength  has  been  calculated,  the  pseudorange 
measurement  noise  error  can  be  established.  The  pseudorange  measurement  noise  error 
values  are  developed  in  a  laboratory  setting  by  testing  the  operating  hardware,  and  vary 
between  GPS  receivers.  The  data  is  compiled  into  a  table  of  corresponding  signal 
strengths,  as  illustrated  in  Figure  3.7.  This  plot  is  expressed  in  received  power  and  not 
C/No.  By  subtracting  No  (-202  dBW)  from  the  signal  strength  values  listed,  the  C/No 
values  can  be  obtained.  This  plot  is  representative  of  an  extremely  high-sensitivity 
receiver.  Similar  C/No  data  is  located  in  [4], 


Figure  3.7:  Pseudorange  Measurement  Noise  Error  Standard  Deviation  Plot  [11] 

The  pseudorange  measurement  noise  error  value  is  given  as  a  standard  deviation, 
which  is  in  turn  multiplied  by  a  MATLAB®  random  number,  providing  the  measurement 
error  that  will  be  added  to  the  pseudorange  value.  The  pseudorange  measurement  noise 
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error  can  be  turned  on  or  off  in  the  simulation.  If  the  error  is  turned  off,  the  pseudorange 
value  will  equal  the  exact  distance  between  the  GPS  and  GEO  satellites,  just  as  in  the 
simple  GPS  measurement  model. 


3.4.4  TWTT  Measurements 

The  Two-Way  Time  Transfer  measurements  in  this  simulation  do  not  include 
Sagnac  error  or  motion  related  errors,  since  they  are  largely  deterministic  and  can  be 
removed.  The  simulation  could  add  the  errors  and  then  remove  them,  but  it  is  a  wasted 
step  that  will  only  increase  computational  cost  and  would  have  no  added  value.  For 
simplicity,  this  simulation  assumes  that  the  propagation  delays  will  cancel  as  in  the  static 
TWTT  case.  The  resulting  TWTT  measurement  equation  is: 

AT  =  y  [TIC  (GEO)  -  TIC  (REF)]  +  oTWTT 
=  ~[2cStGEO  —  2cStREF  J  + 

UTwtt  (3.19) 

—  cStGEO  —  cSt  REE  +  vTWTT 

where 

St GEO  =  GEO  satellite  clock  error 
StREF  =  reference  clock  error 
oTWTT  =  TWTT  measurement  error 
c  =  speed  of  light 

If  TWTT  measurements  are  to  be  used  in  the  simulation,  they  are  given  to  the 
Kalman  filter  along  with  the  pseudoranges  for  measurement  incorporation.  The  TWTT 
measurement  error  can  be  turned  on  and  off,  using  values  of  10,  3,  0.3,  0.03,  and  0.003 
meters. 


55 


3.5  Kalman  Filter 


A  Kalman  filter  was  chosen  over  a  least  squares  batch  filter,  as  was  used  in  the 
research  conducted  in  [8].  The  Kalman  filter  allows  for  the  use  of  new  measurement  data 
as  it  becomes  available  and  easily  allows  for  any  stochastic  processes,  such  as  clock 
errors.  The  least  squares  estimation  algorithm  can  not  use  data  as  it  is  collected, 
requiring  all  information  up  front.  It  also  is  more  difficult  to  model  stochastic  processes 
in  a  least  squares  estimator. 

The  Kalman  filter  has  several  initial  values  that  govern  the  estimation  algorithm. 
The  filter  must  know  the  accuracy  level  of  the  incoming  measurements  and  how  much  to 
trust  in  their  positioning  information,  as  well  as  the  amount  of  process  noise  in  the 
system. 


3.5.1  Initial  State  Values 

The  first  thing  needed  by  the  Kalman  filter  is  an  initial  state,  including  position, 
velocity,  and  clock  error.  These  initial  values  are  gathered  from  the  truth  data,  so  that  the 
filter  will  begin  its  estimation  at  the  same  initial  point  of  the  truth  data.  The  state  vector 
is  given  in  Equation  3.20. 

x  =  [x  Y  Z  X  Y  Z  cStGEO  cStGEO  c5tref  cSiref  ]T  (3.20) 

where 

X,Y,Z  =  GEO  satellite  position  components 

X,Y,Z  =  GEO  satellite  velocity  components 

StGEO  and  SiGEO  =  GEO  satellite  clock  bias  and  clock  drift 

8tref  and  8iref  =  TWTT  reference  clock  bias  and  clock  drift 

c  =  speed  of  light 

If  initial  state  error  is  turned  on,  a  particular  error  standard  deviation  is  multiplied 
by  a  random  number  and  added  to  each  initial  value.  The  initial  state  error  standard 
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deviations  are  listed  in  Table  3.2.  These  were  considered  reasonable  values  for  the  start 
of  a  simulation  run. 

Table  3.2:  Initial  State  Error  Standard  Deviation  Values 


Initial  State  Value 

Standard  Deviation  Value 

Position 

20  m 

Velocity 

0.01  m/s 

Clock  bias 

14  m 

Clock  drift 

20  m/s 

The  dynamics  matrix  F  describes  the  motion  of  the  satellite,  and  is  modeled  using 
the  same  equations  as  the  truth  generation  (see  Equation  3.6).  The  Kalman  filter  F  matrix 
includes  the  clock  terms,  making  it  a  10  x  10  matrix,  as  seen  below. 

x  =  Fx  (3.21) 

o  o  o  iooooo  olr  x  ' 

0  0  00100000  Y 

0  0  00010000  2 

-f  0  0  0000000  x 

0  -4  0  0000000  Y 

r  •  (3.22) 

o  0-40000000  2 

r 

0  0  0  0000100  c5tGEo 

0  0  0  0000000  c5iGEo 

0  0  0  0000001  cStref 

0  0  0  000000  oJL  c5iref  _ 

where 

X,Y,Z  =  GEO  satellite  position  components 

X,Y,Z  =  GEO  satellite  velocity  components 

X,Y,Z  =  GEO  satellite  acceleration  components 

StGEO ,  SiGEO ,  and  StGEO  =  GEO  satellite  clock  bias,  drift,  and  drift  rate 

8tref ,  8iref ,  and  8iref  =  TWTT  reference  clock  bias,  drift,  and  drift  rate 

c  =  speed  of  light 
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Since  the  truth  generation  and  Kalman  filter  use  the  same  dynamics,  they  will 
produce  the  same  results  if  all  errors  are  turned  off.  This  provides  a  method  for 
quantitative  error-checking  of  the  code. 

The  initial  covariance  matrix  P  describes  the  accuracy  of  the  state  vector,  and  will 
be  updated  as  the  filter  iterates.  The  first  six  values  along  the  diagonal  are  for  the  X,  Y,  Z 
position  and  velocity  values.  The  last  four  values  along  the  diagonal  are  for  the  clock 
bias  and  drift  of  the  GEO  satellite  and  TWTT  reference  clocks.  The  P  matrix  uses  the 
values  shown  in  Table  3.2  and  is  described  in  Equation  3.23,  where  subscripts  p  = 
position,  v  =  velocity,  cb  =  clock  bias,  and  cd  =  clock  drift. 

(7p2  0  0  0  0  0  0  0  0  0 

0  cr/  0  0  0  0  0  0  0  0 

00  cr/  0000  0  0  0 

0  0  0  ctv2  0  0  0  0  0  0 

0000  cr2  00  0  0  0 

0  0  0  0  0  crv2  0  0  0  0 

0  0  0  0  0  0  (Jcb2  0  0  0 

0  0  0  0  0  0  0  (jJ  0  0 

0  0  0  0  0  0  0  0  crj  0 

0  0  0  0  0  0  0  0  0  aj 

The  covariance  matrix  Q  describes  the  errors  associated  with  propagating  the 
state  covariance  matrix  P  through  time.  The  Q  matrix  includes  the  process  noise  value  of 
the  GEO  satellite  and  the  related  clock  q  values.  The  process  noise  value  is  equal  to  the 
process  noise  value  used  in  the  truth  generation  and  the  clock  q  values  are  taken  from 
Table  3.1,  depending  upon  the  type  of  clocks  that  are  used.  Equation  3.24  describes  the 
Q  matrix  used  in  the  Kalman  fdter,  and  it  does  not  change  throughout  the  simulation. 
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000000  0  0  0  0 

000000  0  0  0  0 

000000  0  0  0  0 

0  0  0  <7„  0  0  0  0  0  0 

0  0  0  0  cr  0  0  0  0  0 

Q  = 

0  0  0  0  0  <7„  0  0  0  0 

0  0  0  0  0  0  cqfsv)  0  0  0 

0  0  0  0  0  0  0  cq2(sv )  0  0 

000000  0  0  cqfref)  0 

000000  0  0  0  cq2  (ref) 

where 

<yn  =  GEO  satellite  process  noise 
qx  (,s'v)  =  GEO  satellite  clock  bias  process  noise 
q2(sv)  =  GEO  satellite  clock  drift  process  noise 
qx  (ref)  =  reference  clock  bias  process  noise 
q2(ref)  =  reference  clock  drift  process  noise 
c  =  speed  of  light 

3.5.2  Filter  Execution 

Once  the  initial  state  vector,  dynamics  matrix,  and  covariance  matrices  are 
defined,  the  Kalman  filter  is  ready  to  begin.  The  first  step  is  to  calculate  the  state 
transition  matrix  <D  by  using  the  dynamics  matrix  F  and  the  time  interval  At ,  shown  in 
Equation  3.25. 

<J>  =  e FAl  (3.25) 

The  next  step  is  to  calculate  Qd,  the  discrete  time  version  of  the  continuous-form 
Q  matrix. 

At 

Qd  =  j0»(r)Q  <t>T(T)dr  (3.26) 

0 

Next,  the  filter  propagates  the  state  and  covariance,  as  seen  below: 
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x(4)  =  ®(4-4-i)i(Ci) 

P(4 )  =  0(4  -4_1)P(4+_l)or(4  -4-0+Q, 


(3.27) 

(3.28) 


Measurement  information  for  a  Kalman  filter  must  follow  the  equation  listed 

below. 

z(()  =  h[x(4),(]  + v(4)  (3-29) 

where 

z(t  )  =  measurement  vector  at  time  U 

v((.)  =  zero-mean  white  Guassian  vector  process  noise  of  strength  R, 
independent  of  process  noise 

The  GPS  pseudorange  measurements  are  of  the  form: 

p  =  ,/(*’“  -  )'  +  (/"  -  +  (z"  - z„, )'  +  cSt„  - cSt“  +  uPR  (3.30) 

Combined  with  Equation  3.29,  Equation  3.30  implies  the  following: 

W(*"  -o’  +(/“  - y,J  +  <*" -z,J  +cStr« -cSf“  (3.31) 


After  the  h  vector  equations  are  written,  the  observational  partial  derivative 
matrix  H  is  constructed.  The  H  matrix  relates  the  linearized  observations  to  the 
estimated  states,  and  is  expressed  in  Equation  3.32. 


H  = 


H2 


H„ 


where 


(3.32) 


rr  _Shl(±)  _SJ12(X) 

~  ’  **  2  ~ 


Sx 


Sx 


;Hn  = 


Shn{±) 

Sx 


(3.33) 
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Sht(x) 

Sx 


(3.34) 


Sh^ix)  8h\(x)  Sh^ix) 

8xi  5x2  8xm 


The  H  matrix  is  of  size  n  by  m,  where  n  is  the  number  of  collected  measurements 
and  m  is  the  number  of  states  in  x  (in  this  simulation,  m  =  10).  For  example,  if  there  are 
two  pseudorange  measurements  and  one  TWTT  measurement  at  a  given  epoch,  the  H 
matrix  and  measurement  vector  z  will  be: 


where 


e 


x 


e 


y 


e 


z 


0  0  0 


0  0  0  1  0  0  0 

0  0  0  1  0  0  0 

0  0  0  1  0  -1  0 


A 


z  = 


Pi 


TWTT 


Xrec~X 


V(x-  -  )’+(/“-  y,,<  f  +  (z“  -z„cf 

_ ZrecZ/l _ 

-x„f+  (/”'  -yr,J+(z“-z„J2 


TT1  -  T  )2 + </"  -  >2 + ^  y 


(3.35) 


(3.36) 


The  error  term  v(6)  in  the  Kalman  filter  measurement  equation  is  controlled  by 
the  covariance  matrix  R.  The  R  matrix  is  specified  for  a  given  measurement  vector.  This 
means  that  each  pseudorange  measurement  will  have  its  own  values  in  the  R  matrix.  The 
R  matrix  is  typically  based  upon  expected  error  statistics,  which  is  based  upon  knowledge 
of  the  system. 

In  the  example  given  above,  the  corresponding  R  matrix  will  be: 
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(3.37) 


0 

0 

2 

CFTWTT_ 

where 

<rra  =  standard  deviation  value  of  pseudorange  measurement  noise 
<JTWTT  =  standard  deviation  value  of  TWTT  measurement  noise 

It  is  now  time  for  the  filter  to  incorporate  the  measurements,  involving  several 
steps.  The  filter  will  calculate  the  residuals,  compute  the  Kalman  gain,  and  update  the 
state  and  covariance. 

The  residuals  are  the  difference  between  the  predicted  state  and  the  observed 
state.  Residuals  are  calculated  using  the  following  equation: 

r  =  z(0-h(f,)  (3.38) 


R  = 


'  PR\ 

0 

0 


0 

J-2 

'p/12 

0 


The  Kalman  gain  determines  the  weighting  factor  that  is  applied  to  the 
measurement  data.  The  Kalman  gain  depends  upon  the  covariance  and  will  instruct  the 
fdter  to  trust  the  measurement  data  more  or  less,  depending  upon  its  accuracy.  The 
Kalman  gain  calculation  is  listed  below. 

K  =  P(4)Hr[HP(t;)Hr+R]  1  (3.39) 

where 

K  =  Kalman  gain  matrix 

p(4")  =  covariance  matrix  propagated  to  next  epoch  (but  before 
measurement  incorporation) 

H  =  observational  partial  derivative  matrix,  which  relates  the  linearized 
observations  (z)  to  the  estimated  states  (x) 

R  =  measurement  noise  covariance  matrix 

Now  the  Kalman  gain  can  be  used  to  update  the  state  and  covariance,  using  the 
following  equations  (stated  previously  in  Chapter  2): 
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(3.40) 


x(4+ )  =  x(4 )  +  K  (z  -  h  [x(4 ),  tk  ]) 

P(C  )  =  P(4  )  -  KHP(4  )  (3.41) 

Lastly,  the  filter  updates  the  F,  O ,  and  Qd  matrices.  In  the  new  epoch,  the  GEO 
satellite  will  have  a  new  radius,  which  changes  the  F  matrix.  The  new  F  matrix  will  then 
change  the  <I>  and  Qd  matrices  using  the  same  equations  listed  above. 

Before  starting  the  next  iteration,  the  current  filter  information  is  saved  into  a 
history  matrix.  When  the  simulation  is  complete,  the  history  matrix  will  contain  the  state, 
covariance,  and  residuals  values  for  each  epoch. 


3.6  Results  Analysis 

Analyzing  the  simulation  results  involves  comparing  the  Kalman  filter  estimated 
state  with  the  true  state.  The  analysis  depends  on  the  simulation  type,  being  either  a 
single  run  or  a  Monte  Carlo  collection  of  runs. 

The  most  important  result  is  the  three-dimensional  positioning  error,  which  will 
be  expressed  as  Mean  Radial  Spherical  Error  (MRSE).  The  MRSE  is  analogous  to  a 
three-dimensional  Distance  Root  Mean  Square  (DRMS)  value.  For  a  Monte  Carlo 
simulation,  the  MRSE  for  a  particular  epoch  is  calculated  by  using  the  following 
equation: 

IlK  +T,2  +zf) 

MRSE  =  U  -  (3 .42) 
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where 


n  =  number  of  simulation  runs 

x.  =  X'me  -  x(llter  =  difference  of  truth  and  filter  X  position  for  epoch  i 

y.  =  Y"'“e  -  Yifilter  =  difference  of  truth  and  filter  Y  position  for  epoch  i 

z.  =  Z\n,e  -  Zfher  =  difference  of  truth  and  filter  Z  position  for  epoch  i 

The  MRSE  value  can  also  be  calculated  by  using  the  standard  deviation  values 
that  exist  in  the  covariance  matrix  for  a  particular  epoch,  as  seen  below. 

\±{t+y!+t)  - 

MRSE  =  - - - *  Ja-  +  a;  +  a;  (3.43) 

where 

ax,ay,az  =  standard  deviation  values  for  the  X,  Y,  Z  coordinates 

As  the  number  of  simulations  n  increases,  the  two  square  root  values  approach  the 
same  value.  Using  a  Monte  Carlo  simulation,  each  run  will  generate  different  position 
values,  but  the  filter-computed  covariance  values  will  be  the  same  for  every  single  run. 
As  a  result,  the  covariance  standard  deviation  values  from  a  single  run  can  replace  the 
position  values  from  hundreds  of  runs  in  a  Monte  Carlo  simulation. 

To  insure  the  covariance  values  are  accurate  and  the  equality  assumption  is  true, 
the  standard  deviation  values  are  first  compared  with  a  Monte  Carlo  simulation  using  a 
high  number  of  runs.  If  the  MRSE  computed  from  the  Monte  Carlo  simulation  closely 
resembles  the  MRSE  computed  from  the  covariance  values  of  a  single  run,  the  model  is 
validated  and  only  a  single  run  is  necessary  to  calculate  the  MRSE  from  that  point  on. 

After  the  MRSE  is  calculated  for  each  time  epoch,  the  Root  Mean  Square  (RMS) 
is  calculated  for  the  entire  collection  of  epochs,  using  Equation  3.44.  The  final  result  is  a 
single  RMS  value  that  depicts  the  level  of  error  in  the  estimation  filter. 
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RMS  = 


Z 


x. 


1=1 

n 


where 

n  =  number  of  epochs 
Xj  =  MRSE  for  epoch  i 


(3.44) 


3.7  Summary 

This  Chapter  describes  the  five  main  functions  of  the  simuation.  The  five 
functions  include  loading  simulation  parameters,  generating  truth  data,  generating 
measurement  data,  running  a  Kalman  filter,  and  analyzing  the  results.  Chapter  4  will 
discuss  the  results  and  analysis  of  the  simulation. 
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4  Results  and  Analysis 


4.1  Introduction 

This  chapter  provides  the  simulation  results  and  analyses.  First,  preliminary 
results  are  explained  that  will  allow  a  better  understanding  of  the  simulation  data.  The 
baseline  results  are  then  described,  giving  a  reference  point  from  which  to  compare  the 
results  of  each  trade  study.  Five  trade  studies  explore  the  effects  of  using  GPS  ephemeris 
data  from  different  days,  using  a  simple  or  complex  GPS  receiver  model,  using  simple  or 
complex  antenna  gain  patterns,  using  different  combinations  of  clock  types  in  the 
receiver  and  TWTT  reference  clock,  and  varying  the  TWTT  measurement  noise  error. 


4.2  Preliminary  Results 

The  primary  purpose  of  this  research  is  to  explore  the  impact  of  using  TWTT 
measurements  on  the  positioning  accuracy  of  GPS  navigation  in  a  geostationary  orbit.  As 
a  first  step  in  validating  the  simulation  model,  several  Monte  Carlo  simulations  were 
performed  to  analyze  the  statistical  results  and  confirm  the  output  was  reasonable. 

4. 2. 1  Graphical  Results 

The  Monte  Carlo  simulations  consisted  of  100  iterations,  and  sample  raw  data  is 
displayed  in  Figure  4.1.  Each  line  represents  the  difference  between  the  true  position  and 
the  Kalman  filter-predicted  position  over  a  single  simulation  run.  Each  simulation  uses 
the  same  parameters  but  will  utilize  a  different  set  of  random  numbers  produced  by  the 
random  number  generator  in  MATLAB®. 
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100  Monte  Carlo  Runs  -  Y-Direction  Position  Error  (ECEF) 
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Figure  4.1:  Monte  Carlo  Measurements  vs  Time  for  100  Iterations  (Standard  GPS  Receiver 

Sensitivity,  no  TWTT) 


The  Monte  Carlo  simulation  results  are  more  appropriately  displayed  by  using  the 
Monte  Carlo  mean,  Monte  Carlo  mean  plus/minus  the  Monte  Carlo  standard  deviation, 
and  the  covariance  calculated  by  the  Kalman  filter,  as  seen  in  Figure  4.2.  The  blue  line 
represents  the  Monte  Carlo  mean  value  at  each  measurement  epoch.  Since  the  Monte 
Carlo  simulation  is  composed  of  100  iterations,  each  point  on  the  blue  line  is  the  mean  of 
100  values.  The  dotted  black  lines  represent  the  sum  and  difference  of  the  Monte  Carlo 
mean  and  the  standard  deviation  of  the  Monte  Carlo  mean.  The  red  lines  represent  the 
covariance  values  computed  by  the  Kalman  filter.  As  the  number  of  simulation  runs 
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approaches  infinity,  the  plots  of  the  Monte  Carlo  mean  plus/minus  the  standard  deviation 
should  exactly  match  the  plots  of  the  filter-computed  covariance. 


Y-Direction  Position  Error  (ECEF) 


Figure  4.2:  100-Run  Monte  Carlo  Measurement  Mean,  Measurement  Standard  Deviation,  and 

Filter-Computed  Covariance 


The  Monte  Carlo  standard  deviation  results  closely  match  the  filter-computed 
covariance  results,  indicating  that  the  covariance  analysis  will  be  representative  of  the 
Monte  Carlo  simulation.  Monte  Carlo  simulations  that  include  TWTT  measurements 
lead  to  the  same  conclusion  as  the  non- TWTT  case. 

The  original  intent  was  to  prove  that  covariance  analysis  is  sufficient  to  analyze 
results  by  comparing  covariance  analysis  results  to  Monte  Carlo  simulations.  If  single 
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runs  are  proved  adequate,  the  massive  computational  time  of  numerous  Monte  Carlo  runs 
can  be  avoided. 

4.2.2  Root  Mean  Square  Results 

As  explained  in  Section  3.6,  the  simulation  results  are  provided  in  terms  of  Root 
Mean  Square  (RMS)  positioning  error.  The  RMS  positioning  error  allows  the  results  to 
be  condensed  into  a  single  number,  providing  a  quick  and  simple  method  for  comparing 
trade  study  results. 

RMS  values  are  normally  calculated  from  a  collection  of  measurement  data 
gathered  from  numerous  simulations  runs,  such  as  a  Monte  Carlo  simulation.  Another 
method  is  to  use  the  filter-computed  covariance  values  calculated  from  a  single 
simulation  run,  as  opposed  to  the  actual  measurement  data.  As  the  number  of  Monte 
Carlo  simulation  runs  approaches  infinity,  the  RMS  values  calculated  by  using 
measurement  data  should  equal  the  RMS  values  calculated  by  using  covariance  values. 
The  filter-calculated  covariance  values  will  be  the  same  for  each  simulation  run,  so 
Monte  Carlo  simulations  are  not  required  when  using  covariance  values  to  calculate  RMS 
values. 

Two  scenarios  were  selected  to  compare  measurement-calculated  RMS  values 
with  filter  covariance-calculated  RMS  values.  An  ovenized  crystal  clock  was  assigned  to 
the  GEO  satellite  receiver  while  using  two  different  GPS  receiver  sensitivity  settings. 
The  three-dimensional  RMS  results  are  displayed  in  Table  4.1.  These  two  scenarios  use 
GPS  ephemeris  data  from  January  1,  2006  and  the  complex  GPS  receiver  model  with 
complex  antenna  gain  pattern  data.  These  settings  are  the  default  for  all  the  simulations 
unless  otherwise  specified. 
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Table  4.1:  Comparison  of  Monte  Carlo  Simulation  and  Single  Run  Simulation,  no  TWTT 


Standard  Sensitivity 
(32  dB-Hz  cutoff) 

High  Sensitivity 
(12  dB-Hz  cutoff) 

Monte  Carlo  3D  RMS 

83.90 

3.88 

Single  Run  3D  RMS 

76.09 

3.83 

%  Difference 


10.26% 


1.31% 


Plotting  the  3D  position  RMS  values  over  time  provides  another  useful  way  of 
comparing  the  Monte  Carlo  measurement  values  to  the  filter  covariance  values  from  a 
single  run.  If  the  filter  is  properly  estimating  the  position  of  the  GEO  satellite,  the  filter- 
computed  covariance  3D  position  RMS  over  time  curve  on  the  plot  should  match  that  of 
the  Monte  Carlo  simulation.  A  plot  of  the  3D  position  RMS  over  time  is  shown  in  Figure 
4.3. 
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3-D  Position  RMS  Values 


Figure  4.3:  3D  Position  RMS  Over  Time  -  Standard  Sensitivity,  no  TWTT 

The  worst-case  scenario  is  a  standard  sensitivity  receiver  with  no  TWTT 
measurements,  displayed  in  Figure  4.3.  In  the  simulation,  the  standard  sensitivity 
receiver  tracks  between  zero  and  four  GPS  satellites  during  a  24  hour  period  (see 
Appendix  A),  resulting  in  poor  measurement  precision.  The  lack  of  GPS  signals  is 
apparent  when  looking  at  the  3D  RMS  plot,  where  the  measurements  curve  differs  from 
the  fdter-computed  covariance  curve. 

Since  few  measurements  are  available  to  the  fdter,  it  must  rely  heavily  on 
estimating  the  satellite  position.  The  covariance  will  predict  how  the  position,  velocity, 
and  clock  errors  will  propagate,  but  the  lack  of  measurements  will  prevent  the  fdter  from 


72 


making  consistent  corrections  to  its  estimation.  This  is  what  causes  the  disparity  between 
the  Monte  Carlo  results  and  the  covariance  results. 

Since  the  Monte  Carlo  measurements  closely  match  the  filter-computed 
covariance  values  and  carry  the  same  trends,  the  model  can  be  considered  valid.  The 
Monte  Carlo  measurement  results  can  now  be  replaced  with  the  covariance 
measurements  from  a  single  simulation  run,  producing  enormous  savings  in 
computational  time. 

4.3  Baseline  Results 

A  series  of  six  scenarios  were  selected  as  the  baseline  for  comparing  against  the 
trade  study  results.  In  the  baseline  analysis,  an  ovenized  crystal  clock  is  used  for  the 
GEO  satellite  GPS  receiver  and  a  Rubidium  clock  is  used  for  the  TWTT  reference  clock. 
All  six  GPS  receiver  sensitivity  levels  are  used  while  switching  the  TWTT  measurements 
on  and  off,  resulting  in  twelve  data  sets.  The  GPS  receiver  sensitivity  levels  dictate  the 
C/No  cutoff  for  the  receiver’s  ability  to  detect  and  use  a  GPS  signal,  as  well  as  the 
associated  pseudorange  error.  The  sensitivity  levels  range  from  32  dB-Hz  (standard)  to  7 
dB-Hz  (ultra  high)  using  5  dB-Hz  increments. 

4.3.1  Plotting  the  Results 

The  plotted  data  provides  a  glimpse  inside  the  simulation  mechanics  and  helps  to 
show  that  the  model  is  valid.  Certain  trends  are  visible  that  point  to  tangible  reasons  for 
the  nature  of  the  data. 

In  situations  where  the  receiver  has  a  lower  sensitivity,  there  will  be  periods  of 
time  when  no  GPS  satellites  are  visible  and  the  filter  will  need  to  rely  completely  upon  its 
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process  model  estimating  the  state  of  the  satellite.  Over  time,  errors  will  propagate  and 
the  filter’s  prediction  will  diverge  from  the  true  position,  as  will  the  filter’s  confidence  in 
the  accuracy  of  the  predicted  position  values.  When  a  satellite  comes  into  view,  the  filter 
will  once  again  receive  measurement  data  and  can  finally  reduce  the  uncertainty.  In  the 
case  of  clock  bias  error,  these  rapid  corrections  are  plainly  visible  in  Figure  4.4.  The  blue 
line  in  Figure  4.4  represents  the  difference  between  the  filter-predicted  state  and  the  true 
state,  while  the  red  line  represents  the  filter-computed  covariance.  If  the  filter  is  accurate, 
the  measurement  values  should  stay  within  the  boundaries  of  the  covariance. 


GEO  Clock  Bias  Error 


Figure  4.4:  GEO  Satellite  Clock  Bias  Error  Filter  Covariance  for  a  Single  Run  (Standard  Sensitivity 

with  no  TWTT  Measurements) 
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Figure  4.5:  Visible  GPS  Satellites  from  GEO  Satellite  (Standard  Sensitivity) 


In  Figure  4.4,  the  rapid  covariance  growth  is  directly  linked  to  the  periods  where 
zero  GPS  satellites  are  visible.  The  filter-computed  covariance  of  clock  bias  error  grows 
when  there  are  no  available  measurements,  since  the  filter  is  continually  losing 
confidence  in  the  prediction.  When  measurements  become  available,  the  filter  can 
immediately  correct  the  prediction  of  the  clock  state,  shrinking  the  covariance  drastically 
as  the  filter  is  once  again  confident  in  the  predicted  clock  state. 

When  TWTT  measurements  are  included,  the  filter  never  relies  solely  on 
prediction  for  the  clock  state.  The  reference  clock  provides  bias,  drift,  and  drift  rate 
information  to  the  filter,  allowing  it  to  more  accurately  predict  the  velocity  and  position 
states.  Figure  4.6  is  based  on  exactly  the  same  GPS  satellite  visibility  as  Figure  4.4,  but 
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now  includes  TWTT  measurements.  The  drastic  jumps  in  the  covariance  disappear  since 
the  filter  never  experiences  clock  measurement  blackouts. 
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Figure  4.6:  GEO  Satellite  Clock  Bias  Error  Filter  Covariance  for  a  Single  Run  (Standard  Sensitivity 

with  TWTT  Measurements  Included) 

When  using  a  higher  sensitivity  receiver,  the  filter  is  able  to  more  adequately 


correct  the  estimated  state.  Figure  4.7  illustrates  the  clock  bias  error  of  a  medium 


sensitivity  receiver  with  TWTT  measurements  turned  on  and  off.  The  higher  sensitivity 


receiver  is  able  to  track  more  GPS  signals  and  allow  the  filter  to  predict  the  clock  state 


more  accurately  than  the  receiver  depicted  in  Figure  4.6.  When  TWTT  measurements  are 


included,  the  clock  state  prediction  improves  further  and  the  covariance  is  smaller  than 


the  non-TWTT  case. 
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GEO  Clock  Bias  Error 


Figure  4.7:  GEO  Satellite  Clock  Bias  Error  with  TWTT  Turned  Off  (top)  and  On  (bottom)  (Medium 

Sensitivity) 


A  similar  outcome  is  also  apparent  when  looking  at  the  position  error.  For 
example,  the  X-direction  position  error  is  shown  in  Figure  4.8. 
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X-Direction  Position  Error 
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Figure  4.8:  Geo  Satellite  X-Direction  Position  Error  with  TWTT  Turned  Off  (top)  and  On  (bottom) 

(Medium  Sensitivity) 
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Plotting  the  simulation  results  can  explain  behaviors  that  would  naturally  be 
expected  in  the  real-world  system.  By  visually  verifying  these  events  in  the  data,  the 
simulation  model  can  be  trusted  with  greater  confidence. 

4.3.2  Position  Error 

The  compiled  RMS  position  error  data  for  the  baseline  case  is  listed  in  Table  4.2. 
The  X,  Y,  and  Z  position  errors  are  expressed  as  RMS  values  over  time.  For  example,  all 
the  filter  covariance  X  values  over  the  entire  simulation  run  time  are  squared  and 
summed.  This  value  is  divided  by  the  entire  number  of  measurement  epochs  to  calculate 
the  mean  of  the  squared  values.  The  square  root  of  the  mean  value  results  in  the  RMS 
value,  condensing  the  entire  simulation  run  into  a  single  number  for  evaluation.  After  the 
initial  state  is  declared,  the  simulation  runs  for  23  hours  and  45  minutes  with  60-second 
measurement  intervals,  resulting  in  1,426  measurement  values. 

The  three-dimensional  RMS  value,  labeled  3D  in  Table  4.2,  is  calculated  in  a 
slightly  different  manner.  When  using  Monte  Carlo  measurement  data,  the  RMS  values 
are  calculated  to  find  the  three-dimensional  position  error  for  each  measurement  epoch 
(refer  to  Section  3.6).  When  using  filter  covariance  values  to  calculate  a  three- 
dimensional  position  error,  a  root  sum  square  (RSS)  is  used  instead  of  the  RMS  (see 
Equation  3.43).  All  X,  Y,  and  Z  values  are  squared  and  summed  and  then  the  square  root 
is  calculated.  The  result  is  1,426  RSS  values  that  represent  the  three-dimensional 
position  error  for  each  epoch.  In  order  to  reduce  this  data  into  a  single  number  for 
comparison,  the  RMS  value  of  all  the  RSS  values  is  calculated,  providing  the  single 
number  to  describe  the  three-dimensional  position  error. 
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Table  4.2:  Simulation  Baseline  Results  -  RMS  Position  Error  Values  per  Reciever  Sensitivities 


GPS 

Receiver 

Sensitivity 

Cutoff 

(dB-Hz) 

RMS  Position 
Type 

RMS  Position  Error  (m) 

%  Decrease  in 
Position  Error 
by  Using  TWTT 

no  TWTT 

Rb  TWTT 

Standard 

32 

3D 

76.09 

60.32 

20.73% 

X  (Radial) 

66.44 

55.49 

16.48% 

Y  (In-track) 

27.69 

17.60 

36.43% 

Z  (Cross-track) 

24.68 

15.78 

36.06% 

Standard 

Plus 

27 

3D 

31.87 

24.10 

24.36% 

X  (Radial) 

30.19 

22.32 

26.07% 

Y  (In-track) 

8.10 

6.97 

14.02% 

Z  (Cross-track) 

6.18 

5.85 

5.40% 

Medium 

22 

3D 

11.34 

7.71 

32.04% 

X  (Radial) 

10.54 

6.56 

37.76% 

Y  (In-track) 

3.11 

2.97 

4.41% 

Z  (Cross-track) 

2.80 

2.74 

1.97% 

Medium 

Plus 

17 

3D 

6.32 

4.08 

35.43% 

X  (Radial) 

5.87 

3.38 

42.49% 

Y  (In-track) 

1.70 

1.66 

2.38% 

Z  (Cross-track) 

1.60 

1.58 

1.09% 

High 

12 

3D 

3.83 

2.40 

37.43% 

X  (Radial) 

3.56 

1.95 

45.15% 

Y  (In-track) 

1.02 

1.00 

1.65% 

Z  (Cross-track) 

0.97 

0.96 

0.83% 

Ultra  High 

7 

3D 

2.36 

1.48 

37.55% 

X  (Radial) 

2.19 

1.18 

46.14% 

Y  (In-track) 

0.65 

0.64 

1.18% 

Z  (Cross-track) 

0.63 

0.62 

0.68% 

4.4  Trade  Study  1:  Ephemeris  Date 

To  ensure  the  simulation  does  not  depend  upon  the  ephemeris  data  for  a  specific 
day,  nine  additional  days  were  selected  for  comparison.  One  day  was  selected  out  of  each 
year  from  1997  to  2006,  providing  a  comprehensive  evaluation  pool.  Each  simulation 
was  identical,  other  than  the  different  ephemeris  date,  and  used  the  worst-case  scenario  of 
a  standard  sensitivity  receiver  with  no  TWTT  measurements.  The  results  of  the  study  are 
listed  in  Table  4.3. 
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Table  4.3:  Various  Dates  of  Ephemeris  and  Resulting  RMS  Position  Errors 


Ephemeris 

Date 

RMS  Position 
Type 

RMS  Position 
Error  (m) 

%  Difference 

from  l-Jan-06 

1 -Jan-06 

3D 

76.09 

■■■ill 

X  (Radial) 

66.44 

Y  (In-track) 

27.69 

IHllI 

Z  (Cross-track) 

24.68 

■■■ 

27-Mar-05 

3D 

68.39 

10.12% 

X  (Radial) 

61.29 

7.76% 

Y  (In-track) 

22.68 

18.09% 

Z  (Cross-track) 

20.19 

18.20% 

16-May-04 

3D 

70.95 

6.76% 

X  (Radial) 

63.51 

4.41% 

Y  (In-track) 

25.56 

7.67% 

Z  (Cross-track) 

18.61 

24.60% 

6-Apr-03 

3D 

74.29 

2.36% 

X  (Radial) 

64.02 

3.63% 

Y  (In-track) 

28.90 

4.38% 

Z  (Cross-track) 

24.19 

1.98% 

15-Sep-02 

3D 

76.78 

0.91% 

X  (Radial) 

66.42 

0.03% 

Y  (In-track) 

26.67 

3.68% 

Z  (Cross-track) 

27.79 

12.59% 

Ephemeris 

Date 

RMS  Position 
Type 

RMS  Position 
Error  (m) 

%  Difference 

from  1 -Jan-06 

11 -Nov-01 

3D 

79.07 

3.91% 

X  (Radial) 

68.01 

2.36% 

Y  (In-track) 

30.58 

10.44% 

Z  (Cross-track) 

26.28 

6.49% 

6-Aug-00 

3D 

78.18 

2.74% 

X  (Radial) 

68.13 

2.54% 

Y  (In-track) 

30.41 

9.84% 

Z  (Cross-track) 

23.36 

5.37% 

20-Jun-99 

3D 

77.45 

1.78% 

X  (Radial) 

66.36 

0.12% 

Y  (In-track) 

28.74 

3.81% 

Z  (Cross-track) 

27.71 

12.27% 

26-Jul-98 

3D 

80.72 

6.09% 

X  (Radial) 

68.38 

2.92% 

Y  (In-track) 

31.59 

14.08% 

Z  (Cross-track) 

29.03 

17.61% 

10-Aug-97 

3D 

89.50 

17.63% 

X  (Radial) 

71.87 

8.18% 

Y  (In-track) 

40.56 

46.49% 

Z  (Cross-track) 

34.64 

40.34% 

For  additional  insight,  the  mean  value  and  standard  deviation  were  calculated 
across  all  ten  ephemeris  dates,  given  in  Table  4.4. 


Table  4.4:  RMS  Position  Error  Measurement  Mean  and  Standard  Deviation  Across  all  Ten 

Ephemeris  Dates 


RMS  Position 

Mean  (in) 

Standard 

Type 

Deviation  (in) 

3D 

77.14 

5.72 

X  (Radial) 

66.44 

2.97 

Y  (In-track) 

29.34 

4.75 

Z  (Cross -track) 

25.65 

4.58 

The  three-dimensional  RMS  position  error  mean  value  is  approximately  one 
meter  greater  than  the  value  for  the  default  January  1,  2006  ephemeris  date.  This  ensures 
that  the  2006  date  will  be  a  good  representation  of  the  expected  values  for  the  simulation 
at  any  given  time.  The  3D  RMS  standard  deviation  is  less  than  six  meters  which  is  well 
within  acceptable  boundaries  for  instilling  confidence  in  the  simulation  model. 
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Trends  in  the  GPS  satellite  visibility  provide  another  look  at  the  simulation 
model’s  validity.  Ephemeris  data  from  each  of  the  ten  dates  produced  very  similar 
visibility  trends,  as  they  all  vary  from  zero  to  four  with  an  average  of  approximately  two. 
Samples  of  GPS  satellite  visibility  plots  are  located  in  Figure  4.9.  For  the  case  of  a 
standard  sensitivity  receiver  and  no  TWTT  measurements,  the  number  of  visible  satellites 
ranged  from  zero  to  four.  The  only  exception  was  when  five  GPS  satellites  were  visible 
for  a  brief  period  while  using  the  2003  ephemeris  data.  This  was  the  only  time  more  than 
four  satellites  were  visible. 

By  examining  the  RMS  position  error  values  and  visible  GPS  satellite  plots,  it  is 
apparent  that  differing  ephemeris  dates  do  not  have  a  significant  impact  on  the  results  of 
the  simulation.  The  simulation  is  accurate  and  consistent  for  any  ephemeris  data  that  is 
used. 


Number  of  Visible  Satellites  vs.  Time  Number  of  Visible  Satellites  vs.  Time 


Number  of  Visible  Satellites  vs.  Time 


Number  of  Visible  Satellites  vs.  Time 


Figure  4.9:  GPS  Satellite  Visibility  Plots  for  a  GEO  Satellite  Using  Four  Different  Ephemeris  Dates 
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4.5  Trade  Study  2:  GPS  Receiver  Models 

The  simulation  contains  a  simple  and  a  complex  GPS  receiver  models.  The 
simple  reciever  model  is  strictly  concerned  with  satellite  look  angles  and  does  not 
calculate  signal  path  loss,  varying  degrees  of  pseudorange  measurement  noise  error,  or 
varying  gain  patterns.  If  the  GEO  satellite  is  within  a  GPS  satellite’s  main  beam  field  of 
view  and  is  not  obscured  by  the  Earth,  the  receiver  is  assumed  to  be  able  to  gather  a 
pseudorange  measurement.  This  measurement  is  simply  the  vector  distance  between  the 
two  satellites  with  the  addition  of  white,  Gaussian  measurement  noise  with  a  five-meter 
standard  deviation. 

The  complex  receiver  model  calculates  the  received  signal  strength  taking  into 
account  transmit  and  receive  antenna  gain  patterns  and  path  loss,  as  described  in  Section 
3.4.3.  The  received  signal  strength  dictates  whether  or  not  a  pseudorange  measurement  is 
available,  and  if  so,  measurement  noise  error,  as  described  in  Section  3.4.3. 

The  simple  GPS  receiver  model  is  dependent  upon  the  look  angle  of  the  GPS 
satellite.  A  smaller  look  angle  will  limit  the  number  of  satellite  that  the  GEO  satellite  can 
track  simultaneously,  while  a  larger  look  angle  will  allow  the  GEO  satellite  to  track  more 
GPS  satellites.  The  first  part  of  this  trade  study  examines  varying  the  look  angle  of  the 
GPS  satellites  and  the  results  are  seen  in  Table  4.5.  The  42.6°  angle  is  a  common 
representation  of  the  GPS  Block  II/IIA  main  beam,  as  seen  in  [1]  and  others.  The  38°  and 
120°  angles  were  chosen  to  mimic  the  Block  HR  main  beam  and  Block  HR  UHF 
crosslink  beam  as  shown  in  [20],  The  50°  angle  was  chosen  arbitrarily  to  represent  a 
modest  increase  over  the  commonplace  42.6°  Block  II/IIA  main  beam. 
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Table  4.5:  Simple  Model  RMS  Position  Error  Results  with  Pseudorange  Measurement  Noise  Error 

Standard  Deviation  of  5  meters 


GPS  Look 

RMS  Position 

RMS  Position  Error  (m) 

Decrease  in  Position 

An  ale 

Tvpe 

no  T  WIT 

Rh  TWTT 

Error  using  TWTT 

3D 

90.45 

6S.34 

24.45% 

W 

X  (Radial) 

73.93 

61.55 

16.30% 

Y(In-track} 

40.63 

20.65 

49.10% 

Z  {Cross-track} 

32.64 

20.34 

37j69% 

3D 

72.66 

56.84  ! 

21.75% 

42  6’ 

X  (Radial) 

63.93 

52 23 

15.31% 

Y  (In-track) 

27.66 

15.89 

42,56% 

Z  (Cross-track) 

20.67 

15.83 

23.44% 

3D 

51.39 

3920  i 

23.72% 

50s 

X  (Radial) 

47.17 

36,12 

23.43% 

Y  (In-track) 

15.67 

10.56 

32,59% 

Z  (Cross-track) 

13.05 

10.98 

15.85% 

3D 

6.72 

4.34  j 

35.45% 

120= 

X  (Radial) 

6.36 

3.78 

40.65% 

Y  (In-track) 

1.52 

1.49 

2.05% 

Z  (Cross-track) 

1.53 

1.53 

0.41% 

When  comparing  with  the  complex  GPS  receiver  model,  there  are  similarities 
between  the  simple  model’s  42.6°  look  angle  and  the  complex  model’s  standard 
sensitivity  level.  The  120°  look  angle  scenario  is  also  similar  to  the  medium  plus 
sensitivity  level.  These  traits  are  listed  in  Table  4.6. 


Table  4.6:  Comparison  of  Simple  and  Complex  GPS  Receiver  Model  RMS  Position  Error  Values 


GPS  Look 

RMS  Position 

RMS  Po  sitio  n  E  rro  r  (m) 

%  Decrease 

Angle 

Type 

no  TWTT 

Rb  TWTT 

Using  TWTT 

3D 

72.66 

56 S4 

21.78% 

42.6= 

X  {Radial"; 

63  93 

52,23 

IB. 31% 

Y  {In-track) 

27.66 

15.59 

42,56% 

Z  (Cross-track) 

20.67 

15. S3 

23.44% 

3D 

6.72 

4.34 

35.45% 

120= 

X  {Radial) 

6.36 

3.7S 

40.65% 

Y  {In-track) 

1.52 

1.49 

2.05% 

Z  (Cross-track) 

153 

153 

0.41% 

GPS  Receiver 

Cutoff 

RMS  Position 

RMS  Position  Error  (m) 

°/b  Decrease 

Sensitivity 

(dB-Hz) 

Type 

™  T  WTT 

Rb  TWTT 

Using  T’fiTT 

3D 

76.09 

60.32 

20.73% 

Standard 

32 

X  {Radial) 

66.44 

55.49 

16.4S% 

Y  {In-track) 

27.69 

1760 

36.43% 

Z  (Cro  5  s-track) 

24.6S 

15.73 

36,06% 

3D 

6.32 

4.0B 

35.43% 

Medium  Pins 

1  7 

X  {Radial) 

5.B7 

3.3S 

42,49% 

Y  {In-track) 

1.70 

1.66 

2,33% 

Z  (Cross-track) 

1.60 

1JS 

1.09% 

The  approximated  Block  II/IIA  antenna  main  beam  angle  in  the  simple  model 
produces  results  very  similar  to  the  standard  sensitivity  level  in  the  complex  model.  This 
suggests  that  the  simple  model  might  be  sufficient  for  examining  scenarios  that  involve 
nominal  equipment  and  values,  where  it  is  safe  to  assume  all  pseudorange  measurement 
noise  errors  follow  a  standard  deviation  of  five  meters.  Specialized  hardware  with  higher 
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sensitivity  levels  will  most  likely  require  the  complex  model  to  represent  the  system. 
Still,  the  fact  that  the  standard  case  shows  correlation  between  the  simple  and  the 
complex  models  further  proves  that  the  simulation  is  reliable. 

The  second  part  of  this  trade  study  explores  various  pseudorange  measurement 
noise  error  standard  deviation  values  and  their  effect  on  three-dimensional  positioning 
accuracy  for  the  simple  receiver  model.  The  pseudorange  measurement  noise  error 
standard  deviation  is  varied  between  1,  2,  3,  5,  and  7  meters,  the  results  being  listed  in 
Table  4.7. 

Table  4.7:  Simple  Model  RMS  Position  Error  Results  with  Varying  Pseudorange  Measurement  Noise 
Error  Standard  Deviations  and  GPS  Look  Angle  of  42.6° 


Pseudo  range 

RMS  Position 

RMS  Position  Error  (m) 

%  Decrease  in  Position 

Measurement  Error  cf 

Type 

no  TWIT 

Kb  TH  TT 

Error  usin^  TWTT 

3D 

51.20 

3422 

33.17% 

1  m 

X  (Radial) 

44.37 

30.50 

31.24% 

Y  fin-track) 

2101 

10.67 

49.20% 

Z  (Cross-track) 

14.57 

1126 

22.74% 

3D 

59.73 

42.14 

29.44% 

2  m 

X  (Radial) 

52.42 

3S.1S 

27.16% 

Y  (Tn-track) 

23.27 

12.34 

46.97% 

Z  (Cross-track) 

16 6S 

12.90 

22.70% 

3D  ! 

65.46 

4B.26 

26.27% 

3  m 

X  (Radial) 

57.6B 

44.09 

23.55% 

Y  (In-track) 

25.00 

13.6S 

45.26% 

Z  (Cross-track) 

IS. 2  5 

14.05 

23  00% 

3D 

12.66 

56.  B4 

21.78% 

5  m 

X  (Radial) 

63.93 

52.23 

IS. 31% 

Y  (In-track) 

27  -56 

15.  S9 

42,56% 

Z  (Cross-track) 

20.67 

15. B3 

23.44% 

3D  ! 

77.17 

62,37 

19.17% 

7  m 

X  (Radial) 

67.45 

57.25 

15.12% 

/  JTl 

Y  (In-track) 

29.  BS 

17.72 

40.72% 

Z  (Cross-track) 

22.62 

17.28 

23.61% 

There  are  distinct  similarities  between  the  7  meter  standard  deviation  value  in  the 
simple  model  and  the  standard  sensitivity  level  in  the  complex  model.  The  correlation  is 
shown  in  Table  4.8. 
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Table  4.8:  Comparison  of  Simple  and  Complex  GPS  Receiver  Models  Based  on  Pseudorange  Error 


Pseudo  range 

RMS  Position 

RMS  Position  Error  (m) 

Q/6  Decrease 

Measurement  Error  ct 

Type 

no  TWIT 

Rt  rwrr 

Using  TWIT 

3D 

77.17 

62.37 

1917% 

7  m 

X  (Radial) 

67.45 

57.25 

1112% 

Y  (In-track) 

29. S£ 

17.72 

40.72% 

Z  (Cross-track) 

22.62 

17.2S 

23  61% 

GPS 

Cutoff 

RMS  Position 

RMS  Position  Error  (m) 

fyia  Decrease 

Receiver 

(dB-Hz) 

Type 

no  TWIT 

Rb  TWIT 

Using  TWIT 

3D 

76.-09 

60.32 

20.73% 

Standard 

32 

X  (Radial) 

65.44 

3149 

16.4£% 

Y  (In-track) 

27.69 

17.60 

36.43% 

Z  (Cross-track) 

24.6  £ 

117£ 

36.06% 

The  standard  sensitivity  level  in  the  complex  GPS  receiver  model  limits  the 
number  of  GPS  satellites  that  are  visible  to  the  GEO  satellite.  The  few  satellites  that  are 
available  will  most  likely  be  large  look  angles  where  the  antenna  gain  patterns  will 
provide  the  least  amount  of  signal  strength  to  the  receiver.  Low  received  signal  strength 
will  incur  a  higher  pseudorange  measurement  noise  error.  For  the  complex  model,  the 
pseudorange  measurement  noise  error  is  15  meters  at  the  receiver’s  signal  strength  cutoff 
boundary  (i.e.,  the  lowest  signal  power  that  can  provide  a  pseudorange  measurement).  It 
makes  sense  that  a  7-meter  pseudorange  error  value  would  best  represent  the  standard 
sensitivity  level  when  using  the  simple  model. 


4.6  Trade  Study  3:  Gain  Pattern  Models 

As  stated  in  Section  3.4.3,  there  are  two  GPS  antenna  gain  pattern  models 
available  for  the  complex  GPS  receiver  model,  taken  from  two  different  sources.  The 
simple  antenna  model  provides  gain  pattern  data  from  nadir  to  43°,  while  the  complex 
model  provides  data  from  0°  to  180°  off-boresight.  Both  gain  pattern  models  replicate  the 
GPS  Block  II/IIA  transmit  antenna.  The  primary  and  secondary  lobes  are  contained  in 
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the  0°-43°  portion  of  the  transmit  beam  (refer  to  Figure  3.4).  Additional  lobes  exist  at 
higher  angles,  but  the  gain  is  very  small.  Because  the  gains  are  so  small  at  high  angles, 
the  complex  gain  pattern  becomes  more  beneficial  as  the  GPS  receiver  sensitivity 
increases.  Table  4.9  displays  the  results  of  this  trade  study. 


Table  4.9:  Comparison  of  Simple  and  Complex  GPS  Gain  Pattern  Models  and  Resulting  RMS 

Position  Error  Values 


Complex  GPS  Gain 

Simple  GPS  Gain 

Pattern  Model 

Pattern  Model 

°/o  Difference  in  Complex  and 

GPS  Receiver 

Cutoff 

RMS  Position 

RMS  Po  s  itio  n  E  rro  r  (m) 

RMS  P  o  s  itio  n  E  rro  r  (m) 

Simple  Position  Error 

Sensitivity 

(dB-Kz] 

Type 

no  TWTT 

Rb  TWTT 

no  TWTT 

Rb  TWTT 

no  TWTT 

Rb  TWTT 

3D 

76.09 

60.32 

S5.34 

69.54 

12,16% 

15.30% 

Standard 

32 

X  (Radial) 

66.44 

55.49 

72,32 

63.09 

8.85% 

13.69% 

Y  (In-track-) 

27.69 

17.60 

34.84 

20.67 

25.84% 

17.46% 

Z  (Cross-track-) 

24.68 

15. 7S 

25. 98 

20.71 

17.42% 

31.24% 

3D 

31.B7 

24.10 

30.02 

21.57 

-5.80% 

-10.49% 

Standard  Plus 

27 

X  (Radial) 

30.19 

22,32 

28,29 

19.64 

-6.30% 

-12,03% 

Y  (In-track-) 

S.10 

6.97 

S.09 

6.94 

-0.20% 

-0.37% 

Z  (Cross-track] 

6.  IS 

5.S5 

5.95 

5.63 

-3.84% 

-3.73% 

3D 

1134 

7.71 

16.79 

10.98 

48.03% 

42,47% 

Medium 

22 

X  (Radial-) 

10.54 

6.56 

16.04 

9.92 

52.18% 

51.25% 

Y  (In-track-) 

3.11 

2.97 

3. S3 

3.5S 

23.40% 

20.60% 

Z  (Cross-track-) 

2.  SO 

2.74 

3.12 

3.04 

11.68% 

10.82% 

3D 

6.32 

4.0S 

10.82 

6.40 

71.24% 

56.97% 

Medium  Plus 

17 

X  (Radial-) 

5 .87 

3.3S 

1041 

5.74 

77.33% 

70.09% 

Y  (In-track) 

1.70 

166 

2,25 

2,15 

32,70% 

30.06% 

Z  (Cross-track] 

1.60 

1 .58 

l.SS 

1.S4 

17.72% 

16.24% 

3D 

3. S3 

2.40 

7.06 

3. SO 

84.35% 

58.61% 

High 

12 

X  (Radial) 

3.56 

195 

6.  S3 

3.38 

91.71% 

73.05% 

Y  (In-track] 

1.02 

1.00 

136 

1.32 

33.88% 

31.60% 

Z  (Cross-track] 

0.97 

0.96 

1 16 

1  13 

19.44% 

17.69% 

3D 

2,36 

1.48 

4.73 

2.39 

100.30% 

62.08% 

Ultra  High 

X  [Radial] 

2,19 

US 

4.59 

2.12 

110.13% 

79.6m 

Y  (In-track] 

0.65 

0.64 

0.S7 

0.S4 

34.14% 

32.06% 

Z  (Cross-track] 

063 

0.62 

0.75 

0.73 

1991% 

15,10% 

As  mentioned  previously,  the  complex  gain  pattern  is  required  when  the  GPS 
receiver  sensitivity  increases.  For  the  standard  sensitivity,  where  the  small  gains  that 
exist  at  high  look  angles  cannot  be  utilized,  the  simple  gain  pattern  will  most  likely  be 
sufficient  for  accurate  simulation  results. 

When  using  a  “standard  plus”  sensitivity  level,  the  simple  gain  pattern  appears  to 
be  more  accurate  than  the  complex  gain  pattern.  This  is  most  likely  due  to  slight 
differences  in  the  lobe  sizes  and  positions  between  the  two  gain  pattern  models.  The 
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minimum  gain  value  between  the  primary  and  secondary  lobes  is  different  by  a  few 
degrees  between  the  two  models,  which  could  cause  the  accuracy  of  the  simple  gain 
pattern  model  to  be  exaggerated  in  this  particular  situation. 


4.7  Trade  Study  4:  Clock  Type  Selection 

Three  types  of  clocks  are  available  in  the  simulation  model.  An  ovenized  crystal 
clock  represents  an  inexpensive,  commonly-used  clock  that  would  be  onboard  most 
common  satellites.  The  Cesium  and  Rubidium  atomic  clocks  represent  hardware  that 
provides  much  greater  precision  and  accuracy  in  time  measurements.  GPS  time  is 
governed  by  high-quality  atomic  clocks  like  the  Rubidium  clock  in  this  model. 

This  trade  study  examines  the  impact  of  atomic  clocks  on  three-dimensional 
position  error.  Because  pseudorange  measurement  accuracy  is  bound  to  the  quality  of  the 
clock  involved,  it  is  expected  that  higher-accuracy  clocks  will  provide  smaller  position 
errors.  The  results  of  this  trade  study  are  listed  in  Table  4.10. 


Table  4.10:  Comparison  of  Clock  Types  and  Resulting  3D  RMS  Position  Error  Values  (in  meters) 


GEO 

Clock 

TWTT  Ref 
Clock 

Standard 
(32  dB-Hz) 

Standard+ 
(27  dB-Hz) 

GPS  Receive 

Medium 
(22  dB-Hz) 

jr  Sensitivity 

Medium+ 
(17  dB-Hz) 

High 

(12  dB-Hz) 

Ultra  High 
(7  dB-Hz) 

Crystal 

None 

76.09 

31.87 

11.34 

6.318 

3.830 

2.363 

Crystal 

Cesium 

60.33 

24.12 

7.722 

4.103 

2.431 

1.518 

Crystal 

Rubidium 

60.32 

24.10 

7.707 

4.080 

2.396 

1.476 

Cesium 

None 

60.37 

24.16 

7.772 

4.135 

2.447 

1.518 

Cesium 

Cesium 

60.31 

_ 24.10 

7.702 

4.080 

2.400 

1.482 

Cesium 

Rubidium 

60.31 

24.09 

7.697 

4.071 

2.386 

1.464 

Rubidium 

None 

60.36 

24.14 

7.757 

4.112 

2.413 

1.476 

Rubidium 

Cesium 

60.30 

24.09 

7.696 

4.070 

2.385 

1.461 

Rubidium 

Rubidium 

60.30 

24.09 

7.694 

4.067 

2.381 

1.456 
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As  expected,  the  inclusion  of  TWTT  measurements  provides  smaller  position 
errors,  as  does  a  more  accurate  clock.  As  the  sensitivity  of  the  GPS  receiver  increases, 
there  is  a  smaller  percent  decrease  in  position  error.  With  the  high  sensitivity  receiver, 
there  are  enough  quality  measurements  to  drastically  reduce  the  clock  error,  so  the 
addition  of  TWTT  measurements  is  less  significant.  The  percent  decrease  in  3D  RMS 
position  error  per  scenario  is  listed  in  Table  4.1 1. 


Table  4.11:  Percent  Decrease  in  3D  RMS  Position  Error  Compared  to  Scenario  Using  Ovenized 
Crystal  GEO  Satellite  Clock,  no  TWTT,  Standard  Sensitivity  Receiver 


GEO 

Clock 

TWTT  Ref 
Clock 

Standard 
(32  dB-Hz) 

Standard+ 
(27  dB-Hz) 

GPS  Receive 

Medium 
(22  dB-Hz) 

tr  Sensitivity 

Medium+ 
(17  dB-Hz) 

High 

(12  dB-Hz) 

Ultra  High 
(7  dB-Hz) 

Crystal 

None 

20.71% 

35.06% 

36.54% 

35.77% 

Crystal 

Cesium 

24.31% 

31.91% 

Crystal 

Rubidium 

20.73% 

24.36% 

32.04% 

35.43% 

37.43% 

37.55% 

Cesium 

None 

20.66% 

24.19%^ 

31.47% 

34.55% 

36.12% 

35.75% 

Cesium 

Cesium 

20.74% 

24.38%^ 

32.08% 

35.43% 

37.34% 

37.28% 

Cesium 

Rubidium 

20.74% 

24.39% 

32.13% 

35.57% 

37.70% 

38.04% 

Rubidium 

None 

20.68%_ 

24.23% 

31.60% 

34.92% 

37.00% 

37.54% 

Rubidium 

Cesium 

20.75% 

24.40% 

32.14% 

35.58% 

37.74% 

38.16% 

Rubidium 

Rubidium 

20.75% 

24.40% 

32.15% 

35.62% 

37.83% 

38.37% 

It  is  interesting  to  note  that  the  inclusion  of  an  atomic  clock  makes  a  remarkable 
impact  on  the  position  error.  As  long  as  there  is  at  least  one  atomic  clock  involved  in  the 
navigation  process,  the  highest  accuracy  is  achieved.  For  example,  when  a  Rubidium 
TWTT  reference  clock  is  introduced  to  an  ovenized  crystal  clock  on  the  GEO  satellite, 
the  effect  is  the  same  as  having  a  single  Rubidium  clock  on  the  GEO  satellite,  or  a 
Rubidium  clock  in  both  locations.  The  difference  between  these  three  scenarios  is  on  the 
order  of  millimeters  of  three-dimensional  position  error.  This  means  that  a  cheap,  simple, 
durable  clock  can  be  placed  on  a  satellite  in  a  harsh  unstable  environment,  while  the 
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expensive,  complex,  fragile  atomic  clock  can  reside  at  a  TWTT  reference  clock  station  on 
Earth  in  a  controlled  environment. 


4.8  Trade  Study  5:  TWTT  Measurement  Noise  Error 

TWTT  measurements  are  subject  to  measurement  noise  error,  just  as  the 
pseudorange  measurements  are.  It  is  expected  that  as  technology  advances,  TWTT 
measurements  will  become  more  accurate,  lowering  the  amount  of  error.  This  trade 
study  examines  the  impact  of  varying  TWTT  measurement  noise  error,  ranging  between 
standard  deviation  values  of  10,  3,  0.3,  0.03,  and  0.003  meters.  Table  4.12  displays  the 
results  of  varying  the  TWTT  measurement  noise  error  standard  deviation.  Gaps  in  the 
table  exist  because  some  scenarios  were  deemed  of  lesser  importance  and  eliminated  to 
save  computation  time. 


Table  4.12:  Comparison  of  Various  TWTT  Measurement  Noise  Error  Standard  Deviation  Values 
and  Resulting  3D  RMS  Position  Errors  (in  meters) 


GEO 

Clock 

TWTT 

Ref 

Clock 

TWTT 

Measurement 

Error  a  (m) 

Standard 
(32  dB-Hz) 

Standard+ 
(27  dB-Hz) 

GPS  Receive 
Medium 
(22  dB-Hz) 

?r  Sensitivity 

Medium+ 
(17  dB-Hz) 

High 

(12  dB-Hz) 

Ultra  High 
(7  dB-Hz) 

10 

63.582 

8.155 

3.106 

Wmmmmi 

3 

60.911 

24.563 

7.849 

4.224 

2.597 

1.739 

Crystal 

Cesium 

0.3 

60.330 

24.118 

7.722 

4.103 

2.431 

1.518 

0.03 

60.323 

24.111 

7.718 

4.099 

2.424 

1.509 

0.003 

60.323 

7.718 

2.424 

10 

63.575 

8.143 

3.093 

3 

60.899 

24.551 

7.836 

4.202 

2.567 

1.710 

Crystal 

Rubidium 

0.3 

60.317 

24.103 

7.707 

4.080 

2.396 

1.476 

0.03 

60.310 

24.096 

7.703 

4.076 

2.390 

1.466 

0.003 

60.310 

7.703 

2.390 

10 

60.306 

7.700 

2.393 

3 

60.305 

24.090 

7.696 

4.070 

2.385 

1.462 

Rubidium 

Rubidum 

0.3 

60.304 

24.089 

7.694 

4.067 

2.381 

1.456 

0.03 

60.303 

24.088 

7.694 

4.067 

2.380 

1.456 

0.003 

60.303 

7.694 

2.380 
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Varying  the  TWTT  measurement  error  has  a  small  effect  on  the  position  accuracy 
in  the  simulation.  When  the  3D  RMS  position  error  is  24  and  60  meters,  having  a  TWTT 
measurement  accuracy  greater  than  3  meters  will  gain  50-60  centimeters  of  position 
accuracy.  When  the  position  error  is  1-7  meters,  increasing  the  TWTT  measurement 
accuracy  beyond  3  meters  gains  0-30  centimeters. 

Overall,  the  accuracy  of  the  TWTT  measurements  is  not  extremely  critical.  In 
some  cases  it  will  provide  half  a  meter  of  increased  accuracy,  but  it  is  usually  less  than 
that. 

4.9  Additional  TWTT  Study  -  TWTT  Reference  Clock  Locked  to  GPS  Time 

The  investigation  of  the  impact  of  TWTT  measurements  on  GPS  navigation  led  to 
another  scenario  concept  that  involved  a  reference  clock  that  is  synchronized  with  GPS 
time.  If  the  reference  clock  could  exactly  match  GPS  time,  then  the  receiver  clock  error 
could  be  essentially  removed  from  the  pseudorange  measurement  using  the  TWTT 
measurement.  In  this  case,  the  TWTT  reference  clock  errors  are  not  estimated  by  the 
Kalman  filter  (since  they  are  zero),  as  they  were  in  the  previous  cases.  Using  the  TWTT 
with  the  reference  clock  synchronized  to  GPS  time  significantly  reduced  the  position 
error,  as  shown  in  Table  4.13.  The  use  of  TWTT  in  this  mode  resulted  in  60-70% 
improvement  in  positioning  accuracy,  depending  on  the  receiver  sensitivity. 
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Table  4.13:  Comparison  of  3D  Position  RMS  Error  (meters)  When  Using  a  TWTT  Reference  Clock 

that  is  Synchronized  with  GPS  Time 


GEO 

Clock 

TWTT 

Standard 
(32  dB- 
Hz) 

GP 

Standard+ 
(27  dB- 
Hz) 

S  Receivei 

Medium 
(22  dB- 
Hz) 

r  Sensitivity 

Medium+ 
(17  dB- 
Hz) 

r 

High 
(12  dB- 
Hz) 

Ultra 

High 

(7dB- 

Hz) 

Crystal 

None 

76.09 

31.86 

11.34 

6.32 

3.83 

2.36 

Rubidium 

None 

60.36 

24.14 

7.76 

4.11 

2.41 

1.48 

Crystal 

Rubidium 

60.32 

24.10 

7.71 

4.08 

2.40 

1 .48 

Rubidium 

Rubidium 

60.30 

24.09 

7.69 

4.07 

2.38 

1.46 

Crystal 

GPS 

Time 

22.93 

9.01 

4.07 

2.33 

1.44 

0.95 

4.10  Summary 

This  chapter  outlined  the  simulation  results  and  provides  analysis  of  the  data. 
The  baseline  results  provide  the  fundamental  information  used  to  compare  against  each 
trade  study.  The  trade  studies  included  using  ephemeris  data  from  different  days,  using  a 
simple  and  complex  GPS  receiver  model,  using  a  simple  and  complex  GPS  antenna  gain 
pattern,  varying  the  clock  types  found  on  the  GEO  satellite  and  TWTT  reference  station, 
and  adjusting  the  TWTT  measurement  error. 

The  results  indicate  that  changing  the  ephemeris  date  has  no  effect  on  the 
simulation  data,  so  the  model  is  valid  for  any  day.  The  simple  GPS  model  and  GPS  gain 
pattern  proved  useful  for  low  sensitivity  receivers,  but  the  complex  model  and  gain 
pattern  are  required  for  mid-  to  high-sensitivity  receivers.  As  long  as  a  single  high- 
accuracy  atomic  clock  is  available  somewhere  in  the  system,  either  on  the  satellite  or 
supplied  via  TWTT  measurements,  the  position  accuracy  is  the  same.  Having  a  TWTT 
measurement  error  standard  deviation  smaller  than  3  meters  had  little  effect  on  the  three- 
dimensional  position  error  of  the  satellite.  A  value  of  30  centimeters  proved  to  be  as 
beneficial  as  a  value  of  3  centimeters.  Most  importantly,  providing  a  TWTT  reference 
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clock  which  is  slaved  to  GPS  time  offers  radical  position  accuracy  enhancements,  far 
beyond  what  is  achievable  with  current  GPS  receivers  alone. 

Chapter  5  will  provide  summary  conclusions  and  recommendations  for  future 
research  in  this  area. 
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5  Conclusions  and  Recommendations 


5.1  Significance  of  Research 

The  goal  of  this  research  was  to  evaluate  the  benefits  of  Two-Way  Time  Transfer 
measurements  when  used  to  augment  GPS  navigation  at  very  high  altitudes  like 
geostationary  orbit.  At  this  altitude,  the  receiver  is  above  the  GPS  satellite  constellation, 
a  location  where  the  GPS  system  was  not  designed  to  provide  navigation  information. 
The  satellite  geometry  at  GEO  severely  limits  the  number  of  usable  GPS  satellites,  since 
the  Earth  occludes  most  of  the  GPS  signals  from  the  GEO  satellite.  In  order  to  produce 
navigation  information,  a  minimum  of  four  satellites  are  needed  to  solve  for  the  three 
dimensional  position  and  clock  error.  When  less  than  four  satellites  are  available,  the 
clock  error  cannot  be  resolved  and  the  information  provided  is  not  very  useful. 

TWTT  provides  a  method  for  resolving  the  clock  error  when  less  than  four  GPS 
satellites  are  visible.  TWTT  is  a  technique  that  involves  simultaneously  exchanging 
signals  between  two  or  more  clocks  and  is  one  of  the  most  accurate  ways  to  compare 
clocks.  With  TWTT  measurements,  the  clock  error  between  the  GPS  satellite  transmitter 
and  the  GEO  satellite  receiver  can  be  resolved,  allowing  useful  information  to  be  gleaned 
from  any  number  of  visible  satellites,  even  if  it  is  less  than  four. 

This  capability  allows  any  standard  GPS  receiver  to  operate  effectively  on  a  GEO 
satellite  with  reasonable  accuracy.  Accurate  GPS  navigation  in  high-altitude  orbits 
provides  numerous  opportunities,  such  as  automated  station-keeping  in  a  GEO  orbit. 
Also,  by  substituting  automation  and  removing  the  ground-based  ranging  systems,  the 
cost  reduction  incurred  by  reducing  ground  support  is  considerable. 
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To  summarize  the  key  results  into  a  single  figure,  five  scenarios  were  selected:  1) 
Ovenized  crystal  GEO  clock,  no  TWTT,  2)  Rubidium  GEO  clock,  no  TWTT,  3) 
Ovenized  crystal  GEO  clock,  TWTT  to  Rubidium  reference,  4)  Rubidium  GEO  clock, 
TWTT  to  Rubidium  reference,  5)  Ovenized  crystal  GEO  clock,  TWTT  to  clock 
synchronized  to  GPS  time.  Each  of  these  scenarios  were  run.  The  three-dimensional 
Root  Mean  Square  (3D  RMS)  position  error  for  each  simulation  run  is  provided  in  Figure 
5.1. 


Figure  5.1:  Keystone  Plot  -  3D  RMS  Position  Error  vs.  GPS  Receiver  Sensitivity  Levels  and  Clock 

Configurations 


This  is  the  key  plot  that  captures  the  essence  of  this  thesis  research,  as  it  clearly 
describes  the  benefit  of  including  TWTT  measurements  with  GPS  pseuodrange 
measurements  for  navigation.  Introducing  an  atomic  clock  provides  some  positioning 
improvement,  but  the  breakthrough  occurs  when  the  TWTT  reference  clock  is 
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synchronized  with  GPS  time,  offering  a  60-70%  reduction  in  three-dimensional  position 


error. 

The  precise  clock  does  not  necessarily  need  to  reside  on  the  satellite,  but  rather 
can  be  placed  on  Earth  where  it  can  be  better  protected  and  maintained.  Also,  if  the 
precise  clock  fails,  the  on-orbit  satellite  is  not  rendered  inoperable.  Ovenized  crystal 
clocks  are  less  complex  and  more  robust,  making  them  a  better  candidate  for  space  flight, 
further  emphasizing  the  importance  of  keeping  the  precise  atomic  clock  on  the  ground. 

The  idea  of  using  TWTT  measurements  in  GPS  navigation  has  been  explored  in 
the  original  user  environment  on  Earth,  and  provided  exceptional  accuracy 
improvements.  By  now  addressing  TWTT  measurements  used  in  GPS  navigation  in 
high-altitude  orbits,  such  as  GEO,  the  accuracy  improvements  are  even  more  remarkable. 

5.2  Trade  Studies 

A  number  of  additional  trade  studies  were  performed  using  this  simulation. 

Results  are  summarized  in  the  sections  that  follow. 

5. 2. 1  Baseline  Results 

The  baseline  results  are  the  primary  research  results  and  are  compared  with  the 
trade  studies.  These  primary  results  consist  of  six  GPS  receiver  sensitivity  levels  that  are 
used  with  and  without  TWTT  measurements.  The  GEO  satellite  clock  is  an  ovenized 
crystal  clock,  the  TWTT  reference  clock  is  a  Rubidium  clock,  and  the  TWTT 
measurement  error  standard  deviation  is  0.3  meters. 

It  is  apparent  that  increasing  the  sensitivity  level  of  the  GPS  receiver  provides  the 
greatest  improvement  in  measurement  accuracy.  Using  a  single  GPS  receiver  with  a 
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predefined  sensitivity  level,  the  inclusion  of  TWTT  measurements  does  allow  for 
significant  positioning  accuracy  improvements  that  would  otherwise  not  be  available. 
This  improvement  ranges  from  21%  (when  using  a  standard  sensitivity  receiver)  to  38% 
(when  using  an  ultra  high  sensitivity  receiver). 

5.2.2  Trade  Study  1:  Ephemeris  Date 

To  ensure  that  the  simulation  model  does  not  depend  on  ephemeris  data,  ten 
separate  days  were  tested.  The  same  scenario  was  run  for  each  day  and  the  results  were 
examined  to  find  similar  trends  that  would  suggest  correlation. 

The  mean  3D  RMS  position  error  value  is  77.14  meters,  while  the  standard 
deviation  is  5.72  meters.  The  standard  deviation  is  less  than  10%  of  the  mean  value, 
insinuating  that  there  are  no  large  anomalies  when  varying  the  ephemeris  data.  The  GPS 
satellite  visibility  plots  for  each  ephemeris  date  are  very  similar,  each  providing  the  same 
general  level  of  satellite  coverage  to  the  GEO  satellite.  These  findings  indicate  that  the 
output  of  the  simulation  model  will  not  be  corrupted  by  slight  differences  in  ephemeris 
data,  thus  the  model  is  valid  for  any  day. 

5.2.3  Trade  Study  2:  GPS  Receiver  Models 

The  simulation  model  contains  two  types  of  GPS  receiver  models.  The  first  type 
is  a  simple  model  that  only  examines  the  maximum  look  angle  of  the  GPS  satellite  and 
the  cutoff  angle  between  the  GPS  satellite  and  the  Earth  that  obscures  the  GEO  satellite. 
The  simple  model  uses  a  constant  pseudorange  measurement  noise  error  standard 
deviation  value.  The  second  type  is  a  complex  model  that  accounts  for  look  angles,  gain 
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patterns,  signal  strengths,  and  varying  pseudorange  measurement  error  standard  deviation 
values. 

The  first  part  of  this  trade  study  examines  the  effect  of  GPS  satellite  look  angle  on 
the  GEO  satellite  3D  RMS  position  error.  Increasing  the  GPS  satellite  look  angle  has  a 
similar  effect  as  increasing  the  GPS  receiver  sensitivity  level  in  the  baseline  results.  The 
42.6°  look  angle  is  comparable  to  the  standard  sensitivity  level,  while  the  120°  look  angle 
is  comparable  to  the  medium  plus  sensitivity  level. 

The  second  part  of  this  trade  study  examines  the  effect  of  pseudorange 
measurement  noise  error  standard  deviation  on  the  GEO  satellite  3D  RMS  position  error. 
A  standard  deviation  value  of  7  meters  is  comparable  to  the  standard  sensitivity  level  in 
the  baseline  results. 

The  simple  GPS  receiver  model  does  not  have  the  precision  of  the  complex 
model,  but  it  is  still  relatively  accurate.  Using  a  GPS  satellite  look  angle  of  42.6°  and  a 
pseudorange  error  of  7  meters,  the  simple  model  performs  very  similarly  to  the  complex 
model  when  using  an  ovenized  crystal  clock  on  the  GEO  satellite,  with  and  without  a 
Rubidium  TWTT  reference  clock.  This  is  an  additional  step  toward  validating  the 
complex  GPS  receiver  model,  which  provides  the  most  realistic  results. 

5.2.4  Trade  Study  3:  Gain  Pattern  Models 

The  complex  GPS  receiver  model  uses  either  a  simple  or  a  complex  GPS  satellite 
antenna  gain  pattern.  This  trade  study  examines  the  impact  of  various  GPS  satellite 
antenna  gain  patterns  on  the  GEO  satellite  3D  RMS  position  error. 

When  using  a  lower  sensitivity  GPS  receiver,  the  difference  in  3D  RMS  position 
error  values  between  the  simple  and  complex  GPS  satellite  antenna  gain  pattern  models  is 
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the  largest,  but  the  percent  difference  is  the  smallest.  The  3D  RMS  position  error  for  a 
standard  sensitivity  receiver  will  differ  by  9.25  meters,  which  is  12%  of  the  3D  RMS 
position  error  when  using  the  complex  antenna  gain  pattern.  The  reverse  is  true  for  the 
higher  sensitivity  receivers.  The  3D  RMS  position  error  for  an  ultra  high  sensitivity 
receiver  will  differ  by  only  2.37  meters,  which  is  a  100%  increase  changing  from  the 
complex  antenna  gain  pattern  model  to  the  simple  model.  This  trade  study  shows  that  the 
simple  GPS  satellite  antenna  gain  pattern  model  is  accurate  and  correlates  well  with  the 
complex  model,  indicating  the  complex  model  is  reliable  and  accurate  while  also 
providing  more  realistic  results. 

5.2.5  Trade  Study  4:  Clock  Type  Selection 

There  are  three  different  types  of  clock  available  in  the  simulation  model  for  use 
in  the  GEO  satellite  and  the  TWTT  reference  station.  This  trade  study  examines  how 
varying  the  clock  type  in  both  locations  affects  the  GEO  satellite  3D  RMS  position  error. 

It  is  immediately  apparent  that  whenever  an  atomic  clock  is  introduced,  the 
position  results  for  a  particular  GPS  receiver  sensitivity  level  reach  maximum  accuracy. 
If  a  Rubidium  clock  is  placed  on  the  GEO  satellite  and  no  TWTT  measurements  are  used, 
the  end  result  is  the  same  as  having  an  ovenized  crystal  clock  on  the  GEO  satellite  and  a 
Rubidium  clock  as  the  TWTT  reference.  This  trade  study  suggests  that  TWTT  will  allow 
any  satellite  to  benefit  from  an  atomic  clock  without  the  burden  of  carrying  one  on  orbit. 

5.2. 6  Trade  Study  5:  TWTT  Measurement  Noise  Error 

This  trade  study  examines  the  effects  of  altering  the  TWTT  measurement  error 
standard  deviation  value  on  the  GEO  satellite  3D  RMS  position  error.  The  default 
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standard  deviation  value  throughout  the  simulation  runs  is  0.3  meters,  which  is 
representative  of  a  sophisticated  TWTT  setup. 

The  TWTT  measurement  accuracy  does  not  appear  to  have  much  effect  on  the  3D 
RMS  position  error  of  a  GEO  satellite.  When  using  a  standard  sensitivity  GPS  receiver, 
changing  the  TWTT  measurement  error  standard  deviation  from  10  meters  to  3  meters 
gains  approximately  3  meters  of  position  accuracy,  or  4%.  When  using  a  higher 
sensitivity  receiver,  the  3D  RMS  position  accuracy  improvement  is  on  the  order  of 
centimeters.  This  trade  study  indicates  that  a  TWTT  measurement  error  standard 
deviation  value  smaller  than  0.3  meters  gains  relatively  nothing. 

5.2. 7  Additional  TWTT  Study  -  Reference  Clock  Locked  to  GPS  Time 

This  additional  study  addressed  a  TWTT  scenario  where  the  reference  clock  is  set 
to  match  GPS  time.  This  will  allow  the  GPS  receiver  to  maintain  the  exact  time  as  the 
GPS  satellites.  As  a  result,  the  pseudorange  measurements  will  have  all  clock  error 
removed,  allowing  the  calculation  of  a  three-dimensional  position  with  only  three 
satellites,  since  clock  error  no  longer  needs  to  be  solved.  This  configuration  will  allow 
more  than  three  times  the  amount  of  positioning  accuracy  as  compared  to  the  scenario 
using  an  atomic  clock  in  the  TWTT  system. 

5.2.8  Overall  Results 

The  first  three  trade  studies  provide  key  information  that  allows  a  deeper 
understanding  of  the  simulation  model.  By  varying  parameters  and  exploring  the  results, 
the  model  can  be  validated,  ensuring  that  the  output  is  reasonable  and  accurate.  The 
mechanics  of  the  simple  GPS  receiver  model  are  relatively  easy  to  test  and  confirm,  so 
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using  the  simple  model  to  help  validate  the  complex  model  is  a  natural  step  towards  a 
reliable  simulation. 

The  last  two  trade  studies  examine  the  concept  of  aiding  GPS  pseudorange 
positioning  with  TWTT,  and  explore  the  benefits  and  limitations  of  such  a  system.  By 
testing  different  clock  configurations  and  TWTT  measurement  noise  error  values  it  is 
confirmed  that  a  single  atomic  clock  placed  anywhere  in  the  system  can  provide  an 
increase  in  positioning  accuracy,  anywhere  from  21-38%.  By  creating  a  reference  clock 
that  is  synchronized  with  GPS  time,  the  positioning  accuracy  increases  by  70%. 

5.3  Recommendations  for  Future  Research 

5.3.1  Improve  Simulation  Fidelity 

There  are  many  aspects  of  the  simulation  model  that  have  the  potential  for 
enhancement.  This  research  was  focused  on  the  accuracy  improvements  that  TWTT 
measurements  offer,  not  on  complicated  orbital  dynamics,  sophisticated  signal 
transmission  and  reception  mechanics,  or  the  inclusion  of  all  possible  error  sources.  The 
following  paragraphs  describe  portions  of  the  simulation  model  that  could  be  improved  to 
include  many  of  the  items  discussed  above. 

The  simulation  could  be  enhanced  to  model  each  GPS  satellite  appropriately  by 
using  specific  information  on  the  individual  Block  of  GPS  satellite  used  by  each  PRN. 
Currently,  the  simulation  uses  the  same  transmit  antenna  gain  pattern  and  transmit  power 
for  all  GPS  satellites.  Obviously,  every  single  GPS  satellite  in  orbit  is  not  a  Block  II/IIA 
variant,  so  the  simulation  model  should  more  accurately  represent  all  of  the  current 
satellite  variants  in  orbit.  Some  of  the  newer  GPS  satellites  provide  higher  transmit 
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power,  which  could  provide  greater  three-dimensional  positioning  accuracy  than  this 
simulation  indicates. 

Similarly,  the  simulation  should  model  the  actual  GPS  receiver  that  is  built  into 
the  GEO  satellite.  The  simulation  currently  uses  a  low-cost,  off-the-shelf  GPS  receiver 
patch  antenna  that  provides  4dB  gain  or  less.  If  the  GEO  satellite  is  using  a  higher- 
quality  antenna  that  provides  higher  receiver  gain,  the  results  would  also  provide  a 
greater  positioning  accuracy  than  indicated  in  this  research. 

A  sophisticated  propagation  model  (not  just  simple  two-body  dynamics)  for  the 
GEO  satellite  orbit  truth  propagation  and  fdter  propagation  should  be  developed.  The 
simple  propagator  does  not  account  for  perturbations  in  the  Earth’s  gravity  field,  the 
gravity  of  the  sun  and  moon,  solar  wind  pressure,  atmospheric  drag,  or  any  other  high- 
order  phenomenon.  These  effects  could  alter  the  results  of  the  simulation,  and  should  be 
modeled  for  complete  accuracy. 

5.3.2  Investigate  a  Non-standard  GPS  Receiver  Antenna  on  the  GEO  Satellite 

This  simulation  assumes  that  the  GPS  receiver  points  directly  to  the  center  of  the 
Earth  and  has  maximum  gain  at  a  zero  degree  look  angle  (boresight).  From  GEO,  the 
Earth  blocks  most  of  the  GPS  satellites  from  view,  so  all  visible  GPS  satellites  will  be  at 
an  angle  from  the  GEO  satellite.  It  is  intuitive  that  the  GPS  receiver  antenna  should  not 
have  maximum  gain  that  points  toward  the  center  of  the  earth  if  the  received  signals  will 
come  from  an  angle.  Instead,  the  GEO  satellite  GPS  receiver  antenna  should  have 
maximum  gain  at  an  angle  that  will  coincide  with  the  direction  of  the  visible  GPS 
satellites. 
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A  GPS  satellite  transmit  antenna  uses  a  gain  pattern  with  maximum  gain  at 
approximately  10°,  providing  the  maximum  signal  power  to  the  limb  of  the  Earth  where 
atmospheric  loss  will  be  the  greatest  for  a  ground  user.  It  is  recommended  that  the  GEO 
satellite  GPS  receiver  have  a  similar  gain  pattern,  allowing  the  maximum  signal  power 
reception  to  be  in  line  with  the  visible  GPS  satellites.  By  boosting  received  signal  power 
in  this  manner,  the  pseudorange  measurement  noise  errors  will  be  reduced,  providing 
more  accurate  navigation  information. 
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Number  of  Visible  Satellites 


Appendix  A  -  GPS  Satellite  Visibility  Plots 


Number  of  Visible  Satellites  vs.  Time 


Figure  A.l:  GPS  Satellite  Visibility  from  GEO  for  Standard  Sensitivity  Reciever  (32  dB-Hz  cutoff) 


Number  of  Visible  Satellites  vs.  Time 


Figure  A.2:  GPS  Satellite  Visibility  from  GEO  for  Standard  Plus  Sensitivity  Reciever  (27  dB-Hz 

cutoff) 
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Number  of  Visible  Satellites  »  Number  of  Visible  Satellites 


Number  of  Visible  Satellites  vs.  Time 


A.3:  GPS  Satellite  Visibility  from  GEO  for  Medium  Sensitivity  Reciever  (22  dB-Hz  cutoff) 


Number  of  Visible  Satellites  vs.  Time 


Figure  A.4:  GPS  Satellite  Visibility  from  GEO  for  Standard  Plus  Sensitivity  Reciever  (17  dB-Hz 

cutoff) 
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Number  of  Visible  Satellites  vs.  Time 


Figure  A.5:  GPS  Satellite  Visibility  from  GEO  for  High  Sensitivity  Reciever  (12  dB-Hz  cutoff) 


Number  of  Visible  Satellites  vs.  Time 


Figure  A.6:  GPS  Satellite  Visibility  from  GEO  for  Ultra  High  Sensitivity  Reciever  (7  dB-Hz  cutoff) 
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